REAL-TIME PRODUCTION AND SETUP SCHEDULING OF DETERMINISTIC AND STOCHASTIC MANUFACTURING SYSTEMS By MOHSEN EL HAFSI A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 1995 780 n n c #•"" This dissertation is dedicated to my parents, Bou AH and Najet, my loving wife Essia, and my beautiful son, Mohamed Amine. ACKNOWLEDGMENTS I would like to take this opportunity to express my deepest thanks and sincere appreciation to Dr. Boghos Sivazlian, chairman of my supervisory committee, for all his invaluable advice and guidance at both the professional and the social levels. I wish to acknowledge and express my appreciation to Dr. Sherman Bai, who served as cochairman of my supervisory committee and immediate supervisor of my dissertation, for all his invaluable advice, guidance and encouragement throughout this research. His guidance and instruction helped the shape and form of this research, as well as, my personal and career development. I would also like to thank the other members of my committee, Dr. Jack Elzinga, Dr. Chung Lee and Dr. Mark Yang, for their comments and suggestions concerning this research. Thanks are also extended to the Department of Industrial and Systems Engineering at the University of Florida for providing financial support. Finally, I wish to express the most heartfelt gratitude to my dearest wife Essia for her constant support and encouragement throughout my studies. And, to my parents Bou Ali and Najet for all the sacrifices they endured and continue to endure so that I succeed and excel in life. m TABLE OF CONTENTS ACKNOWLEDGMENTS iii ABSTRACT ^ CHAPTERS 1 INTRODUCTION 1 1.1 Motivation 2 1.2 Production and Setup Planning Models 2 1.2.1 The economic order quantity model 3 1.2.2 Extension of the EOQ model: The Wagner-Whitin model 4 1.2.3 Extensions of the Wagner-Whitin model 5 1.2.4 The economic lot scheduling problem 7 1.3 Real-Time Setup Scheduling with Feedback 9 1.3.1 Feedback control scheduling framework 9 1.3.2 Corridor policies n 1.3.3 Distributed real-time setup scheduling policies 11 1.3.4 Other policies and results 12 1.4 Assumptions .13 1.5 Dissertation Outline 14 2 MATHEMATICAL FORMULATION AND SOLUTION APPROACH DETERMINISTIC SYSTEM 15 2.1 Introduction 15 2.2 Problem Formulation 17 2.3 Solution Approach 20 3 OPTIMAL STEADY STATE SOLUTION DETERMINISTIC SYSTEM 27 3.1 Introduction 27 3.2 The Two-Part-Type Case 29 3.2.1 Problem formulation 31 3.2.2 Optimal solution 35 3.2.3 Numerical examples 41 IV 3.2.4 Special cases 43 3.3 The Multi-Part-Type Case 50 3.3.1 Problem formulation 50 3.3.2 Algorithm 52 3.3.3 Numerical example with backlog allowed 55 3.4 Summary 57 OPTIMAL TRANSIENT SOLUTION OF A TWO-PART-TYPE MANUFACTURING SYSTEM WITH SETUP TIMES AND NO SETUP COSTS DETERMINISTIC SYSTEM 58 4.1 Introduction 58 4.2 Optimal Location of the Cyclic Schedule 58 4.3 Partition of the x-Space into Two Mutually Exclusive Regions 60 4.4 Optimal Transient Solution in Region 5R° 63 4.4.1 Initial surplus levels in Region G 65 4.4.2 Initial surplus levels in Region Gl 72 4.4.3 Initial surplus levels in Region G2 73 4.5 Optimal Transient Solution in Region 9t u (Algorithmic Approach) 76 4.6 Optimal Transient Solution in Region 9t u (Analytical Approach) 89 4.6.1 Preliminary discussion and results 89 4.6.2 Equations of Boundaries (Bl and <B2 93 4.6.3 Corridor windows 97 4.6.4 The complete analytical solution 102 4.6.5 Special cases.. 105 4.6.5.1 No backlog for Part Type 2 105 4.6.5.2 No backlog for either part type 107 4.6.5.3 No inventory for Part Type 2 108 4.6.5.4 No inventory for either part type 110 4.7 Summary 112 OPTIMAL TRANSIENT SOLUTION OF A TWO-PART-TYPE MANUFACTURING SYSTEM WITH SETUP TIMES AND SETUP COSTS DETERMINISTIC SYSTEM 113 5.1 Introduction 113 5.2 Optimal Location of the Extracted Cyclic Schedule 113 5.3 Partition of the x-Space into Two Mutually Exclusive Regions 116 5.4 Optimal Transient Solution in Region SR° 117 5.5 Optimal Transient Solution in Region 9t u 122 5.5.1 Equations of Boundaries <M and <B2 124 5.5.2 Corridor windows 124 5.5.3 Numerical example 127 5.6 Case of Zero Setup Times and Nonzero Setup Costs 128 5.6.1 Optimal location of the extracted cyclic schedule 128 5.6.2 Optimal transient solution 131 5.7 Summary 134 6 OPTIMAL AND NEAR OPTIMAL CONTROL OF A TWO-PART-TYPE STOCHASTIC MANUFACTURING SYSTEM WITH NONRESUMABLE SETUPS 136 6.1 Introduction 136 6.2 Mathematical Formulation 136 6.2.1 System dynamics 137 6.2.2 Control constraints 140 6.2.3 Penalty Function 140 6.2.4 The stochastic control problem 141 6.2.5 Optimality necessary and sufficient conditions 142 6.3 Numerical Approach 143 6.3.1 Markov chain approximation method 143 6.3.2 Computational algorithm: The value iteration method 146 6.3.3 Implementation considerations 147 6.4 Optimal Policy 148 6.4.1 Optimal setup policy when the machine is up 149 6.4.2 Optimal setup policy right after repair 150 6.4.3 Optimal production policy 153 6.4.4 Optimal cost function 153 6.4.5 Policy structure 153 6.5 Demand Feasibility 157 6.6 Near Optimal Real-Time Policy for Sufficiently Reliable Systems 160 6.6.1 Construction of the hedging corridor policy 160 6.6.2 Performance of the hedging corridor policy 161 6.7 Summary 167 7 SUMMARY AND SUGGESTIONS FOR FUTURE RESEARCH 168 7.1 Summary 168 7.2 Suggestions for Future Work 170 LIST OF REFERENCES 172 BIOGRAPHICAL SKETCH 176 VI Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy REAL-TIME PRODUCTION AND SETUP SCHEDULING OF DETERMINISTIC AND STOCHASTIC MANUFACTURING SYSTEMS By Mohsen El Hafsi August, 1995 Chairman: Boghos D. Sivazlian Cochairman: Sherman X. Bai Major Department: Industrial and Systems Engineering Many manufacturing systems incur a significant setup time and/ or setup cost whenever the production resource is set up for a new product after it has been used for another. Setups require a downtime during which there is no production; in turn, this might result in decreased production. Cost may be the result of material consumption due to testing (scraps), energy consumption, man hours lost, and tolerance adjustment of the machine for the next part. This situation leads to the following decision making process: when to stop this multi-product system from manufacturing the current product or part, and when to start manufacturing the next product. Basically, the problem is to select a batch size for each product to be manufactured and to sequence the products on the machine. The purpose of this research is to provide real-time optimal or near- optimal feedback policies for determining production rates and setup times for deterministic and stochastic systems (i.e., machine subject to failures and repairs), where the demand of each product is assumed to occur at a constant Vll rate. The most critical issue in real-time scheduling and planning of manufacturing systems is the computational effort required to achieve good performance. To reduce the complexity of the system, most of the quantities of interest are modeled as continuous. We deal with continuous parameters, such as production and demand rates and cumulative production surplus (inventory or backlog) . In order to gain insight and develop solution techniques, we first deal with a two-product system. That is, the real-time production and setup scheduling of a one-machine two-product manufacturing system. The machine can be either perfectly reliable (deterministic system) or unreliable (stochastic system) . The deterministic system is formulated as a feedback control problem, where the state of the system is measured by the amount of cumulative production surplus (inventory if positive and backlog if negative). The performance measure or objective is to minimize the total setup, inventory and backlog costs over a finite planning horizon. The optimal steady state and transient solutions are obtained separately. The stochastic system is formulated as a stochastic control problem. The objective is to minimize the expected total discounted cost of setup, inventory and backlog over an infinite planning horizon. To obtain the optimal policy, a numerical approach is adopted. Based on the results of the deterministic model, a simple heuristic policy, implementable in real-time, is proposed. Simulation results show that the heuristic policy performs very well when applied to manufacturing systems involving low failure rates. vm CHAPTER 1 INTRODUCTION 1.1 Motivation Many manufacturing systems can produce more than one kind of product. In most, a significant cost and/or time is incurred whenever the manufacturing system is set up to produce a new product after it has been used for another. Setups require a down time during which there is no production, which in turn, results in decreased production. Cost may be the result of material consumption, due to testing and tolerance adjustment of the machine for the next part. The reader is referred to the book by Gershwin (1993) for good examples of manufacturing systems involving setup time and/or cost. There are situations where it could be more economical to purchase one machine that is capable of processing many products or parts than to purchase many dedicated machines. This situation leads to the question of how to decide when to stop this multiple-product machine from producing the current product or part, and which product to start making next. Basically, the issue is to select both a sequence in which products will be produced and a batch size for each product to be produced. If setups are performed too often, high change over costs are incurred; and, if too much time is spent on setups, there is not enough time to produce the required amount of the product to meet the demand. On the other hand, if setups are not performed often enough, larger inventory levels of the part being produced may result while shortages may occur for the parts not being produced. The major goal of this research is to provide real-time optimal or near- optimal feedback policies for deciding on production rates and setup times, for both stochastic (subject to failures and repairs) and deterministic systems, in which one is only concerned with operations and setups. In order to gain insight and develop solution techniques, we will be dealing only with small scale systems. However, although small scale, the problems under consideration are still very challenging and practically meaningful. The most critical issue in real-time scheduling and planning of manufacturing systems is the computational effort required to achieve good performance. To reduce the complexity of the system, most of the quantities of interest will be modeled as continuous. Basically, we will be dealing with continuous parameters, such as production and demand rates and production surplus. Once the solution is obtained, it is then easy to convert the production rate into discrete loading times at the machine. In the remainder of this introductory chapter, we will review some models of setup planning without feedback then give an overview of the literature on real-time setup scheduling with feedback. 1 .2 Production and Setup Planning Models There has been an enormous amount of literature on the subject of setups, lot sizes, and related issues. We do not survey all of it; we mostly review the part of the literature that is related to the kind of problems we are dealing with in this research proposal. The models described in this chapter are chosen because of either historical or economical importance. Furthermore, these models form a good starting point for the real-time scheduling models presented in this proposal. In section 1.2.1 we present one of the oldest and simplest lot sizing models, the Economic Order Quantity model. In section 1.2.2, we present an elegant extension of this model due to Wagner and Whitin (the Wagner- Whitin model). In section 1.2.3, the Wagner- Whitin (1958) formulation is extended for the case of multiple parts under capacity restrictions (Newson's model, 1975a, 1975b). In section 1.2.4, the Economic Lot Scheduling Problem, which is a more detailed model, is presented. 1.2.1 The economic order quantity model The simplest model of lot sizing is the Economic Order Quantity or EOQ model. The EOQ model does not describe a true setup manufacturing system because it deals with only a single item or product. It describes a system that acquires items rather than produces them. The items are replenished at discrete times rather than continuously due to ordering costs rather than setup costs. There are many variants of the EOQ formula. A simple version is given by 2KD <? = ,— (1-1) where Q is the economic order quantity, if is a fixed ordering cost incurred each time an order is placed, D is the constant demand rate of the items and h represents the dollar cost per unit of holding an item in inventory. Although this is a very simplified model, the EOQ formula describes a feedback control policy. First, given the quantities K, D, and h, we compute Q. Using Q, the inventory level is monitored as follows: whenever it reaches 0, order Q. Here, the inventory level represents the state of the system. Hence, this is a feedback, or real-time, control policy since the ordering action is based on a measurement of the system state. Other variations of the EOQ model have been studied. These include the effect of lead time, the effect of random demand and the effect of shortage, to cite a few. 1.2.2 Extension of the EOQ model: The Wagner-Whitin model Wagner and Whitin (1958) extend the EOQ model by dropping the assumption of a steady-state demand rate. They divide the time into periods of fixed length. The demand for items is known during each period. The model assumes a single product. Their objective is to minimize the total inventory holding costs. They treat material as a continuous quantity. The dynamics of the system is described by the following difference equation: i-l i-1 where /, is the inventory entering period t (state variable), x t is the amount of material procured or produced during period t (decision or control variable), and d t is the amount demanded during period t. The problem is formulated as a dynamic programming problem with value function f t satisfying /,(/,)= ngn[Wi +S(x t )s t +f t JI t +x t -d l )l 1=1, ..., N; f N Vn ) - mi? [i,J„ +S(x N )s N ] x N 2V I N +X s md N where i t is the interest charge per unit of inventory carried forward to period r+1, s t is the ordering or setup cost, N is the number of periods in the planning horizon and S(x t ) is a binary function given by f if x = 0; S(x) = { [1 if x>0. Notice that this model assumes that the resource has to be set up if an item is manufactured during time period t regardless of whether the same kind of item was manufactured during t - 1 . The optimal solution of this problem suggests that new production occurs only when the inventory becomes zero. The amount produced each time is exactly enough to meet the demand until the next time the inventory goes to zero (this makes sense because the system is deterministic and the demand is known in advance). Because the Wagner-Whitin model assumes perfect knowledge of the demand for future periods (deterministic system), it is appropriate for production planning rather than real-time scheduling. 1.2.3 Extensions of the Wagner-Whitin model A significant extension of the Wagner-Whitin formulation is due to Newson (1975a, 1975b), who explicitly represents multiple parts in the formulation as well as capacity constraints. The objective is to minimize the sum total of setup, production and inventory holding costs. Because of the capacity constraint this problem is called the Capacitated Lot Size Problem (CLSP). As in the Wagner-Whitin formulation Newson divides the planning horizon into T periods of fixed length. He treats material as a continuous quantity. The problem is formulated as a Linear Prograrnming problem. The decision variables are I it , the amount of inventory of Product i at the end of period t, and x H , the amount of Product i made in period r. There are two constraints. The inventory dynamics are given by lit-i+x* = h+<k> ,i'=l, -, N, t = l, ..., N; I i0 > 0, i= 1, ..., N are given where N is the number of products to be produced and d lt is the demand of Product i during period t. The capacity constraint is expressed as ZK**»)+i£*]*n,, fc = i, •••> K, t-1, ..., r i-l where rj* is the capacity lost due to a setup for Product i on resource k, r£ is the per unit capacity lost due to the production of Product i on resource k, R^ is the level of resource k available in period t, K is the total number of resources and S(x it ) has the same meaning as in the Wagner- Whitin formulation. Additional constraints include the nonnegativity of the decision variables I it and x it (i-l, ..., N, t = 1, ..., T). The objective function is given by f = f d f d s i S(x a ) + v i x it +h i I tt . i-l t-1 where s ( is the setup cost for Product i, v i} the per unit production cost of Product i and h { is the inventory holding cost per unit of Product i per period. Like the Wagner-Whitin model, Newson's model is appropriate for production planning rather than real-time scheduling. Although this is a highly simplified formulation, there is no explicit analytic solution to this problem. To solve this problem Newson and others provide heuristics and numerical techniques. Other extensions of the Wagner-Whitin formulation include Mixed Integer Programming formulations. These problems are known to be hard to solve (i.e., their computational time grows exponentially with the size of the problem) . 1.2.4 The economic lot scheduling problem The Economic Lot Scheduling Problem (ELSP) is usually stated as follows: There is one or more identical machines on which m products are to be produced. The demand rate for products is assumed constant, no backlog or backorders are allowed and the system is in steady state. The problem can be viewed as one of finding lot sizes that meet the demand and minimize the average setup and inventory holding costs per unit time. Elmaghraby (1978) gives a thorough review of the models used for the ELSP. Recent research in the field includes the work of Goyal (1984), Dobson (1987), and Carreno (1990). This section is based on Dobson's (1987) formulation. The problem is one of deciding on a cycle of length T, a production sequence F 1 , ..., F n (the sequence may contain repetitions, i.e., n>m) to be repeated indefinitely, production runs *,, ..., t n and idle times iu„ ..., w n so as to meet the demand (no backlog is allowed) and minimize the sum total of inventory holding cost and setup cost. The data for the problem are for each Product i (i =1, ..., m), u, denotes the constant production rate, d { denotes the constant demand rate, h { denotes the inventory holding cost, A, denotes the setup cost and s, denotes the setup time. Dobson's model as well as many others assumes sequence-independent setups. That is, the setup time depends only on the part to be produced next (if it is different from the previously produced) . The decision variables are the position F', in a given sequence where part i is produced; the cycle length T, the length of the production sequence n (since each product is allowed to be produced more than once in a production sequence, n>m); the duration of the i th production run t, ; and the idle time w i during the i production run. Let J,. = {j\F J = i} be the set of times in the cycle, that Product i is produced. Let L k be the set of positions in a given sequence, from k (when F k is produced), up to the next position in the sequence that F k is produced again. The definition of L k assumes that the sequence F 1 , ..., F n is circular. Let 3 be the set of all feasible finite sequences of parts. That is, 3 = {(F 1 , ..., F n )\m < n < oo, \j.\ > 1, i = 1 ,..., m} . Using this notation the ELSP is formulated as follows: ^inf ™m]M)-h F Au Fi -d Fi )^(t pi f + A Fi u.20 J" 1 \ pi subject to X",.* FJ =4T,i=l, ..., m; (1.2) j*j. ±(t FJ+ s FJ+ w Fj ) = T. (1.4) Constraints (1.2) state that the total production over a given sequence must meet the demand for each product. Constraints (1.3) state that the inventory level of Product i reaches zero just as production of Product i is about to begin. This is usually known as the Zero Switch Rule (ZSR) in the literature. The third constraint states that the total amount of production, setup and idle times over a cycle must be the length of the cycle. To solve this problem, Dobson develops a heuristic to determine a production sequence by first obtaining approximate production frequencies. These frequencies are rounded to powers of two. Dobson then uses a bin- packing heuristic to obtain a sequence. Additional cost penalties are incurred when the sequence is determined and the exact quantities and tuning of the batches are computed. Dobson's formulation differs from that of Newson (presented in Section 1.2.3) in that setup behavior is modeled more precisely. It is more restrictive in that the demand rate is required to be constant. Hence, Dobson's model is more appropriate for shorter planning horizon compared to that of Newson. 1 .3 Real-Time Setup Scheduling with Feedback All the models presented in the previous sections, focus on the setup scheduling problem in isolation and under steady-state conditions. They do not address nonsteady-state situations caused by disruptions such as machine failures and scheduled maintenance in a manufacturing system. Under non- steady-state conditions the problem cannot be reduced to one of lot sizing since the lot sizes will not remain constant over time. The purpose of this research is to investigate the setup scheduling problem in a feedback control framework, which could respond to various discrete events in a manufacturing system. 1.3.1 Feedback control scheduling framework Feedback control models for production scheduling appeared in the early '80s, with the pioneering work of Kimemia and Gershwin (1983). They 10 developed a feedback scheduling framework and suggested that the emphasis should be on the flow rate of parts rather than on individual jobs or operations. This framework assumes manufacturing systems that are operated to meet a multiple-part demand rate as closely as possible. It focuses attention on random failures and repairs of machines affecting the production capacity. Based on the works by Rishel (1975) and Olsder and Suri (1980), Kimemia and Gershwin derived a feedback solution to the problem of dispatching parts to machines in a failure-prone Flexible Manufacturing System or FMS (an FMS is a manufacturing system that can produce many different part types with little time and expense lost in changeover from one product to another). Their approach was to separate the relatively long term issues (machine failures) from the short term decisions (part loading). They showed that the optimal control has a special structure, described by a hedging point strategy, in which a positive production surplus (or inventory) of parts is maintained during times when excess capacity is available to hedge against future capacity shortages brought by machine failures. The state of the system is measured by the production surplus of all parts, which is defined as the difference between the cumulative production and the cumulative demand of a part during a certain interval of time. Kimemia and Gershwin's approach is often referred to as flow control in the literature. Recently, Gershwin (1993) proposed a more general framework suggesting a hierarchical approach for the planning and scheduling of manufacturing systems. He groups events from the least frequent (at the top level of hierarchy) to the most frequent ones (at the bottom level). If the hierarchy levels are well separated (frequency wise), each level can be formulated as a flow control problem. 11 1.3.2 Corridor policies Using the formalism of Kimemia and Gershwin (1983), Sharifnia et al. (1991) studied the setup scheduling problem of a multiple part machine. They proposed Corridor Policies to time the setup switching epochs. The corridor policy does this by constructing a corridor around the production surplus trajectory in the production surplus space. The corridor is made up of a set of planes (in n-dimensional space) that are each labeled with the index of a part; when the trajectory hits a plane, the setup is changed so that the machine can produce the part corresponding to the plane just hit by the trajectory. Sharifnia et al. showed that for a round robin sequence (in which the parts are produced only once per sequence and in a cyclic fashion) and when the corridor satisfies certain conditions, the production surplus vector converges to a cyclic schedule which they refer to as the Limit Cycle. Following the ideas of Sharifnia et al., Srivatsan and Gershwin (1990) developed methods for choosing the parameters of the corridors when the setup frequencies are not all the same. 1.3.3 Distributed real-time setup scheduling policies Perkins and Kumar (1989) and Kumar and Seidman (1990) studied distributed deterministic systems. They proposed a class of setup scheduling policies implementable in real-time in a distributed fashion. These policies make setup change based on buffer levels upstream of a machine. To define the proposed policies, consider a system with a single machine that produces N parts. Type i parts arrive at the rate d. at a buffer upstream of the machine. Type i parts need r, time units to be processed on the machine. Then, the policies are defined as follows: 12 Clear the Largest Buffer (CLB): Do not change setup until the buffer of the parts being processed is cleared; then switch to the buffer with the largest fa, (buffer level) . Clear the Largest Work (CLW): Do not change setup until the buffer of the parts being processed is cleared; then switch to the buffer with the largest r,.fa; (work load) . Clear a Fraction (CAF): This is a class of policies which generalizes CLB and CLW. Do not change setup until the buffer of the parts being processed is cleared; then switch to the buffer i for which fa, > Sf . bj , For some s . Perkins and Kumar and Kumar and Seidman showed that under this class of policies the system remains stable. That is, the system can operate with finite buffer levels. Lou et al. (1991) used the same policies to study a deterministic single- machine-multiple part system. Their main result is a characterization of the dynamics induced by CLW policies and tighter bounds for the buffers (compared to the ones provided by Perkins and Kumar, 1989). Lou et al. (1992) extended their results to the case when the machine is unreliable and showed that the total Work-In-Process (WIP) is a recurrent stochastic process. 1.3.4 Other policies and results Connolly (1992) proposed a heuristic for the two-part one-machine setup system, based on known results from the non-setup system. Her approach is based on a local optimization that maximizes the progress toward a target cyclic schedule. Her model considers backlog only. 13 Bai and Burhanpurwala (1993) studied a stochastic two-part, two- machine, two-buffer manufacturing system with dynamic setup changes. They developed a four-level hierarchical real-time production control algorithm. Gallego (1993) studied the ELSP problem in the case of a machine subject to disruptions of small magnitude. He shows that the optimal policy selects the production lot sizes as a linear function of the current inventory levels. Hu and Caramanis (1995) studied a three-part setup system. First, they solved the problem numerically to deduce the structural properties of the optimal policies. Based on the numerical results, they proposed near-optimal policies. The latter are based on primary and secondary setup change surfaces (in n-dimensional backlog/ surplus space). When encountered, primary switching surfaces (which are of dimension n-1) indicate that the current setup must be changed. Secondary switching surfaces (which are of lower dimension) indicate which setup to change to. Hu and Caramanis also proposed a procedure for designing primary and secondary switching surfaces. 1 .4 Assumptions Throughout this proposal, we assume that the following assumptions hold: • The system has infinite supply of raw material • The system has infinite space available for the storage of finished parts; • Many copies of each product or part are made; • The parts are grouped into lots or batches of specified sizes; • The products or parts are not required to be made in a single batch; • We can decide when exactly to perform a setup change; 14 • The system produces perfect parts or products (no rework); • All non-fulfilled orders or tardy jobs are backlogged. 1.5 Dissertation Outline In Chapters 2, 3, 4, and 5, we study a manufacturing system consisting of a single perfectly reliable machine producing two part types. In Chapter 2, we give the general mathematical formulation of the problem and its solution approach. In Chapter 3, we study the system at the steady state. The optimal steady state solution of the two-part-type system is derived analytically. The formulation is then extended to the multi-part-type system and the optimal solution is obtained numerically by means of a modified classical nonlinear programming algorithm. In Chapter 4, we study the transient behavior of the two-part-type system with negligible setup costs. In Chapter 5, we extend the results of Chapter 4 to study the transient behavior of the more general two- part-system including setup times and setup costs. In Chapter 6, we study a system consisting of a single machine subject to random failures and repairs. The setup times are random and nonresumable. That is, after each repair completion, the machine must be set up for a part type. A stochastic optimal control formulation of the problem is given. The problem is then solved numerically. Based on the numerical solution and results from the deterministic case, a heuristic, implementable in real-time, is proposed for the setup switching policy. Simulation results for testing the heuristic performance are provided. In Chapter 7, we summarize our results and indicate possible future directions of this research work. CHAPTER 2 MATHEMATICAL FORMULATION AND SOLUTION APPROACH DETERMINISTIC SYSTEM 2.1 Introduction The goal of the Production and Setup Scheduling Problem (PSSP) is to determine the optimal production rates and setup epochs of several products on a single machine. The latter is assumed to have controllable production rates (i.e., the production rate can be any value between zero and a certain maximum value). Each product has known constant demand rate and processing time. When production is switched from one product to the next a constant setup time as well as a fixed setup cost are incurred. Backlog is allowed and the system is not necessarily in steady state. The objective is to control the production rate of each product as well as to control the setup change epochs so as to minimize the total setup, inventory and backlog costs over a finite or infinite planning horizon. The problem reduces to the famous Economic Lot Scheduling Problem (ELSP), when the planning horizon is infinite, the system is in steadystate, the machine has fixed production rates, no backlog is allowed, and the objective is to determine lot sizes that minimize the average setup and inventory holding costs per unit time. The PSSP is concerned with the study of the system during its transient period (i.e., before stabilizing at the steady state) as well as during its steady state period, throughout the planning horizon. For the deterministic system, we consider a finite planning horizon and assume that the steady state is finite time achievable. This model may have many applications in practice. For 15 16 instance, consider the following situation. We would like to plan the production of several products, periodically (say quarterly), on a bottleneck machine. The demand rate for each period of time is forecasted at the beginning of that period and should be satisfied during the planning period. The machine has enough capacity to meet the forecasted demand for each product and is fast enough to bring the inventory/ backlog of all products to a steady state (which constitutes the most economical way of operating the system) . Now, at the end of the period (planning horizon), the inventory and backlog of the products are checked and (possibly) different forecasted demand rates are to be used in the next period. It is clear that the output (inventory/backlog) of the previous period constitutes the input for the current period, thus new initial conditions as well as new parameters (different demand rate) are to be used. As can be seen, at the beginning of each planning period we are faced with the problem of optimally driving the system to its steady state production cycle. Many other factors can change the conditions of the system and might drive it away from its steady state. Thus, the importance of the transient optimal solution. The PSSP can be formulated as a feedback control problem. The control must respond to certain initial disruptions so as to minimize a certain criterion. This kind of formulation is usually classified under production flow control models. In this chapter, we mathematically formulate and give the solution approach for a deterministic one-machine two-part-type system within a feedback control framework. 17 2.2 Problem Formulation Consider a single-machine manufacturing system producing two distinct part types (or products); each has a constant demand rate d\ [i = 1,2). When production is switched from Part Type j to Part Type i (/*i), a given constant setup time S i and setup cost fc. (i = 1,2) are incurred. Our formulation follows the general framework introduced by Kimemia and Gershwin (1983), where the production flow is modeled as continuous rather than discrete. Let x t (t) be the production surplus of Part Type i (i = 1,2) at time t. A positive value of x t (t) represents inventory while a negative value represents backlog. Let u { (t) be the controlled production rate of the machine producing Type i parts at time t. Let oft) = (c^(t^tir 3 (t),0j |9 (t),0' u (t)) be the setup state vector of the machine at time t. Where, o;.(t), cr s (f) (j#i, i = 1,2, j = 1,2) are right continuous binary functions of t, such that cr. (f ) = 1 when the machine is ready to produce Type i parts and cr(t) = otherwise. er s (f) = l when the machine is undergoing a setup change from Part Type./ to Part Type i and ov (t) = otherwise. Let s(t) be a nonnegative right continuous function of t which takes on the value <^ at the beginning of each setup change to Part Type i (i = 1,2) and decreases with time. s(t) indicates whether a setup is completed or not. We assume that initially the machine is not set up for either part type. System dynamics and constraints The dynamics of the system can be described by ^ = u,.(f) - d,.,i = l,2; (2.1) at 0<u i (t)<[/ i .o;(t),i = l,2 (2.2) where U { is the maximum production rate of the machine when it is producing Type i parts. 18 The setup states of the machine obey the following set of constraints: Oi(t) + tr 2 (t) + <r ia (t) + o- w (f) = 1; (2.3) if o;(r) = 1 and a t [t) = 0, then s(t) = S. and a 9 (t) - 1; (2.4) 'ifs(f )>0 andoj(r) = l, thens(t)=-l ond<r f (t)»l; (2.5) if s[t~) = and ff s (f ) = 1, then ct^ (t) = and s(f) = and oj(t) = 1 (2.6) for i = 1,2, J = 1,2, i # J, where s(t) denotes the time derivative of s(t). The above equations, (2. 3) -(2. 6), merely state that if <y.(t) = l, we can either continue producing Part Type i, or decide to switch production to Part Type /. In the latter case, we must spend exactly Sj amount of time setting up the machine for Part Type j. That is, <j^ (t) = 1 for exactly 8, amount of time (0< s(t)< Sj). After the setup change, tr.(t) = l and the machine is ready to produce Part Typ e J- Penalty function For mathematical convenience, we assume that setup costs are incurred at a constant rate % l =k l /$ l (i = l,2) dollars per unit time during a setup change. Hence, at the end of a setup change to Part Type i, we would have a total setup cost of k { . The instantaneous cost which penalizes the production for being ahead of (i.e., x i > 0) or being behind (i.e., x i < 0) the demand is given by h(x) = ^(c; x ;(t) + c- x :(t)) where c* and c~ are the per unit instantaneous inventory holding and backlog costs respectively, and x*(t) = max{x i (t),0} and x~(t) = max{-x i (t),0}. The total instantaneous cost is then given by Sf^^-hW+^^^^W-S^UXW+^XW+^MKW). (2-7) 19 State variable The state variable of the system is given by x(t ) = (x 1 (t) } x 2 (t)). The state space is given by x e SH 2 . Control variables The control vector u(t) = (UjftJ^ft)) is the vector of production rates of the two part parts. The control vector a{t) = (<J 1 [t),a 2 [t),cr 12 (t),a 21 (t)) is the vector of setup states of the system. We denote by (<r,u) the complete control vector. Production control space The production control space represents the set of all feasible production rates. When the setup state is a{t) at time t, it is given by n(otf)) = {u(t) | < u,(r) < U r a t {t), i = 1,2} . (2.8) Hence, for each setup state we have a different production control space. These are given by Q(0,0) = {u(f)| u,(t) = 0, i = l,2}; £1(1,0) = {u(t) | < u,(t) < U iy i^it) = 0}; (2.9) Q(0,1) = {u(t) | u,(t) = 0, < ^(t) < U 2 }. Usually the production control space is referred to as the capacity set. Setup control space The setup control space is the set of all possible setup vectors o{t) = {cr l (t),cr 2 (t),cr 12 (t),a 21 {t)) satisfying constraints (2.3)-(2.6). Let 3> denote this set. Admissible control policies We denote by E(®,Q.) the set of feasible controls. The set of admissible control policies A is the set of all mappings ju, y: W 2 -* E(n,®), which satisfy /j.(x) = (cr,u) and which are piecewise continuously differentiable. These 20 admissible control policies are feedback controls that specify the control actions (setup and production rate of the machine) to be taken, given the state of the system. Objective function The objective is to determine an optimal control policy jx* e A, corresponding to a setup control a* = {o\,cf 2 ,cf 12 ,a 21 ) and a production flow rate control u* = [ii[,tC,), that minimizes for each initial state x(t) the following cost functional: J,{x(%t) = f g(x(s),a{s))ds (2.10) where, the miriimization is over all functions //(x(r)) = (ofr),u(r)) such that x(r), o[t), and u(t) satisfy the system dynamics (2.1) and (o[t),u{t)) eH(0,f2) for t < t<,t f , where t f is assumed to be sufficiently large. 2.3 Solution Approach The optimal solution to the above problem can be obtained in two parts by dividing the planning horizon (\t,t f \) into two distinct periods. A transient period ([t,t a ]) and a steady state period ([t s ,r r ]). The steady state period corresponds to the case where the state of the system (inventory/ backlog or production surplus levels) has already reached a cyclic schedule, where the produced lots for each part type are of constant size over time. The transient period corresponds to the period where the state of the system still has not reached the cyclic schedule. Formal definitions of the steady state and transient solutions for our problem are given as follows: 21 Definition 2.1: When the production surplus trajectory follows a cyclic schedule in x-space (or production surplus space), we say that the system has reached a steady state. Definition 2.2: An optimal steady state solution corresponds to a cyclic trajectory in x-space minimizing the total average setup, inventory and backlog costs per time unit. Definition 2.3: A transient solution is defined as a trajectory in x-space emanating from an initial point and reaching the steady state in a finite amount of time. Definition 2.4: An optimal transient solution is a transient solution that minimizes the cost of reaching the steady state. That is, among all transient trajectories that lead to the cyclic schedule, it is the one that gives the minimum total cost (setup, inventory and backlog costs). In this section, we state and prove a theorem that will allow us to reduce the set of feasible production rates from an infinite set to a finite one with only three possible production rates. Since we distinguish between the transient and the steady state periods, the cost functional given by (2.10) can be written as the sum of two components: One is due to the transient period, and the other is due to the steady state period. Let t s be the time instant the system reaches the steady state. The total cost can then be written as follows: J,(x(t),t) = £ g(x(s), o{s))ds + £' g(x(s), o{s))ds ; = Jl(x(t),t) + (t f -t s )J*{x(t s ),t s ). (2.11) We refer to J T (x{t),t) as the transient cost, 22 Jl(x(t),t) = £ g(x(s)Ms))ds; (2.12) And J^(x(t s ),t s ) as the average steady state cost, </>(U*J = r-^—yl\{x(s),ais))ds. (2.13) Throughout this dissertation, we assume that the following condition holds. Al: We assume that (t f -t), the planning horizon, is long enough so that the system reaches the steady state and stays there for a long period of time. The following Lemma is based on Assumption Al : Lemma 2.1: Let t s be the time the system reaches the steady state, and Let The the length of the cyclic schedule. Then, we have lim f' g(x(s), a{s))ds = -[g(x(s), a\s))ds; J s And by assumption Al, we have 1 fr J^(x(t s ),t s ) = -lg(x(s)Ms))ds. Proof: Let t f -t s = nT+t , where n is a nonnegative integer. Then, — f g(x(s), o{s))ds = — 1— £ + " r+ ' g(x(s), o{s))ds; — f g(x(s), a{s))ds + —^— f +t ° g(x(s), a{s))ds + 1 Jo nT + 1 Jt ' h •. nT , . giving lim — — \' f g(x(s),o(s))ds = lim — ^— f sr(x(s),o(s))ds + — — P +< ° g(x(s), a(s))ds 23 Since [ g(x(s), o{s))ds and |"° +t °0(;x:(s),a{s))ds are finite it follows that Now by assumption Al , we have ^ (*(*,)>*.) = , lim -(' g(x(s),o{s))ds = — f g(x{s),o{s))ds. ■ The following theorem will be very important for the solution of the deterministic system. Theorem 2.1: The optimal production rate vector u*(s) = (u^sjjU^s)), (t <,s<t f ), belongs to the finite set of vectors fi* = {(0,0), (t/,,0), K,0), [0,U 2 ), (0,dj)}, where u'(s] = if the machine is not producing or undergoing a setup change. The machine produces at the demand rate (Uf(s) = d^) only when x,.(s) = (i = l,2). Proof: First, we prove the theorem for the transient case. Then, we show that a similar argument can be used for the steady state case. The proof is based on the Hamilton- Jacobi-Bellman (HJB) equation. Throughout the proof, we assume that the optimal cost functional is differentiable in x and t. In fact, we will show later in the dissertation that the optimal state trajectory is continuous and piecewise linear. Hence, the optimal cost will not depend explicitly on t and will be the sum of quadratics in x (since the cost rate is linear in x) and therefore differentiable in x. Let J' (x,t) = min \ f g(x(s),o[s))ds, where T f is finite. 24 Transient case: In this case, we let T f = t s in the expression of J' T to obtain the optimal transient cost component. The HJB equation (see Gershwin (1993) for a formal derivation) is given by: dJ' Tj (x,t) = mm dt r >» . f dJ'{x,t) t dJ' T (x,t) } It is clear that, when the machine is undergoing a setup change to a Part Type, there is no decision to make and (u^,u^) is forced to be equal to (0,0). Now, assume that we know the optimal setup state of the machine. Let <r= (1,0,0,0) be this setup state. That is, the machine can produce Part Type 1. In this case, the HJB equation can be rewritten as follows: &T,[x,t) . f , _, SJ' T (x t t) cJUx,t) --5E--JR{**"»+-^h-41+-^h-4jJ. Notice that at each time instant t, if we knew J* (x,t), we would solve a linear programming problem for which u, and u, are the decision variables, dfrj^q and 8J' T Jdx 2 are the cost coefficients and £2(1,0) is the constraints set. Q(1,0) = {(1^,1^)10 <Uj <U lf u^ =0}is bounded and convex. We know that the solution of the above linear programming problem is always at a vertex of the constraint set «(1,0). That is, K,^*) is either equal to (0,0) (if 8J\ Idx^ >0) or equal to (I/„0) (if dJ'^jdx^ < 0). Furthermore, the solution is unique if the cost coefficient MfJ&q is nonzero. In the case dJ\ s Jdx x =Q, the solution is not unique anymore since any (u,*,i^) will not affect the objective function of the linear programming problem at time instant t. However, to keep the cost coefficient dJ^Jdx^ equal to zero at time instant t + St, we should produce Part Type 1 at the demand rate ^ so as to minimize the rate of increase of the cost 25 function J* y In this case (u^,u 2 ) is equal to (d,,0). A similar argument is used when the optimal setup state is cr= (0,1,0,0). Steady-state case: In this case, let t s = t and T f = t f in the expression of J' T above. Let J* = min J s Jx,t) = nun r -vf g(x(s), a{s))ds . Here, we have an average cost formulation. The HJB equation (see Kushner and Dupuis (1992) for a formal derivation) in this case is given by where Wx,t)= lim Jl -(t, -f)J*. As in the previous case, for each time instant t, if we knew V[x,tj, we would solve a linear programming problem which decision variables are the production rates i^ and t^ and which cost coefficients are dV/cb^and dVjdx^ respectively. Hence, using a similar argument as for the transient case and using Lemma 1 , the result follows immediately. ■ The general approach for the deterministic system consists of determining the steady state optimal solution first (which is in principal easier to get than the transient one), which corresponds to an optimal cyclic schedule in x-space. Based on Lemma 2.1, the steady state optimal solution consists of solving a special economic lot scheduling problem (ELSP) where backlog is allowed and the production rates are controllable. To obtain the transient optimal solution, we partition the x-space into mutually exclusive regions. For initial production surplus vectors x(t ) in each region, we determine an optimal 26 trajectory emanating from x(t) and reaching the cyclic schedule (which becomes a target set in this case) with minimum total cost and in finite amount of time. Each region will be characterized by a family of optimal trajectories leading to the cyclic schedule. CHAPTER 3 OPTIMAL STEADY STATE SOLUTION DETERMINISTIC SYSTEM 3.1 Introduction As mentioned in Chapter 2, the steady problem corresponds to a special Economic Lot Scheduling Problem (ELSP) where, in addition, the production rates are controllable and backlog (negative production surplus) is allowed. In the following, we give a brief background of the ELSP due to its importance as a research problem. The ELSP deals with the scheduling of the production of several products on one or more identical machines. Each product has a known constant demand rate and a fixed production rate when it is being produced. When production is switched from one product to the next a sequence- independent constant setup time as well as a fixed setup cost are incurred. The time horizon is infinite, the system is in steady state and no backlog is allowed. The objective is to determine lot sizes that minimize the average setup and inventory holding costs per unit time. In Gallego (1989), the ELSP is extended to allow backlog. In this chapter, we consider the case where nonsatisfied demand is totally backlogged. Although there has been a large amount of research work on the ELSP, an optimal solution approach has not been proposed yet. Rather, good (some times excellent) heuristics have been suggested. A comprehensive review of the ELSP through 1976 is given in Elmaghraby (1978). Recent work on the ELSP includes the work of Boctor (1982), Hsu (1983), Maxwell and Singh (1983), 27 28 Axsater (1984), Goyal (1984), Roundy (1985), Dobson (1987), Gallego (1989), Jones and Inman (1989), Lee and Surya (1989), Carreno (1990), and Zipkin (1991). Most of the aforementioned work emphasized the feasibility of cyclic schedules. A cyclic schedule can be one of two types: Common Cycle schedule or Basic Period schedule. The Common Cycle schedule is also known as the Rotation Cycle schedule, where production is cycled through the products every T units of time. It is well known that this type of schedule is always feasible. Basically, for the Common Cycle schedule, the sequencing problem is eliminated. In the Basic Period schedule, each product's cycle time is an integer multiple of a basic cycle time. Since the production rates are constant, all lots of each product are of equal size. Usually, the Basic Period schedule gives a lower cost than the Common Cycle schedule. The main problem of the Basic Period schedule is its feasibility. In fact, for this kind of approach, Hsu (1983) has shown that even the problem of finding a feasible schedule is an NP-hard problem. To overcome this feasibility problem, Dobson (1987) suggested a new formulation that allows time varying production runs and which includes setup times explicitly in the problem formulation. The main advantage of this approach is that it always provides a feasible schedule. Recently, researchers have dropped the assumption of fixed production rates at the machine and have improved the ELSP model using controllable production rates. The production rate of each product can be chosen to be any value less than a maximum rate. Recent work along this new direction includes the work of Buzacott and Ozkarahan (1983), where they studied the case of a two-product system. First, they categorize the products according to their dollar value of usage, then they showed that only the product with the lower dollar value of usage is produced at maximum. Arizono et al. (1989) studied the effects of controllable production rates on inventory systems. They showed 29 that a controllable production rate inventory system is more efficient than one with fixed production rates. Silver (1990) studied a multi-product system under the Common Cycle schedule assumption. He showed that at most one product slows down its production rate. The optimal production rates and cycle time length were obtained numerically. In all of the aforementioned work, the production rates were decided at the beginning and supposed to be fixed during a product's production run. Moon et al. (1991) generalized Silver's model by considering a system with controllable production rates during the production runs of the products. They obtained the optimal production rates numerically and showed that savings almost twice as large as those reported in the literature can be obtained. The chapter is organized as follows. In Section 3.2, we present the notation and state the assumptions of the model. In Section 3.3, we study the two-product problem. In Section 3.4, we extend the formulation to the multi- product problem and present an algorithm to obtain the optimal solution numerically. We conclude the chapter with Section 3.5. 3.2 The Two-Part-Tvpe Case In Chapter 2, it has been shown that the steady state period corresponds to the case where the state of the system (production surplus levels) has already reached a cyclic schedule, where the produced lots for each part type are of constant size over time. If the machine were perfectly flexible (i.e., with zero setup change times and costs), it would be optimal to produce both parts simultaneously at the demand rates (theoretically, it is optimal to instantaneously switch production back and forth between the two products, since no cost and time are incurred 30 when we switch) . Thus, keeping the production surplus at the zero level. In the case of significant setup times and/or costs, it is not possible to produce both part types at the same time. Thus, at the steady state the production surplus vector must follow a cyclic schedule that repeats itself over time. Based on Theorem 2.1 defined in Chapter 2, it is not difficult to see that a feasible cyclic schedule has the general structure shown in Figure 3.1, where x, represents the inventory /backlog axis of Part Type i (i = 1,2). Now, assume that we start the cycle at Point PI, we then progress toward Point P2 by producing Part Type 1 at the demand rate along segment [P1,P2], where x 2 decreases until we reach Point P2. At this point, we increase the production rate to the maximum and continue producing Part Type 1. x, increases while x 2 decreases until we reach Point P3, where we switch production to Part Type 2. During the setup of the machine for Part Type 2, both inventory levels decrease. Once the machine is ready to produce Part Type 2 (Point P4), the production begins with the maximum allowable rate in order to eliminate backlog as soon as possible (since after the setup, we endup with a backlog for Part Type 2). Once Part TyP e 2 backlog is completely eliminated, we decrease the production rate to the demand rate so that the system moves along [P5,P6\. When we reach P6, the production rate of Part Type 2 is increased to the maximum and a certain inventory is built to hedge against future shortages brought about setups and production of Part Type 1 . At Point P7, we setup the machine for Part Type 1. At Point P8, we produce Part Type 1 at maximum production rate to eliminate backlog as soon as possible until we reach Point PI, where we start the cycle all over again. The steady state solution is completely characterized when the optimal location and shape (in x-space) of the cyclic schedule is known. Using Lemma 2.1 of Chapter 2, it can be seen that determining the cyclic schedule is 31 equivalent to solving a special two-product ELSP. Graphically, the optimization problem can be seen as one of locating and determining the shape of the cyclic schedule in the x-space so as to minimize the average setup, inventory holding and backlog costs per unit time. In the following, we formulate the problem mathematically and then derive the optimal solution in closed form. Notice that we do not consider idle times. The reason is as follows: the purpose of idling the machine (i.e., stopping the machine completely) during the production of a product is usually to stretch the cycle time and delay setup costs which may be high in some cases. In our case, this task is accomplished by producing the products at the demand rate which eliminates inventory and backlog costs for a product during its production at the demand rate along with delaying setup costs. In the case of fixed production rates, during idle times there is a certain inventory or backlog that is present for which a cost is incurred. Therefore idling cannot be optimal in our case. 3.2.1 Problem formulation Notation: For Product i (i = 1,2) t t : time spent producing at maximum rate r, : time spent producing at the demand rate S, : maximum inventory level s { : maximum backlog level y i = cfcT /(c* + c~): cost factor p i = d i /U i : utilization factor of the machine by Product t T = yZ (t, + r, : + S t ): duration of the cyclic schedule 5 = 5 l + S 2 : total setup time during T 32 Figure 3.1: General Structure of the Cyclic Schedule K = it, + fc, : total setup cost during T p = p 1 +p 2 : total utilization factor of the machine a^{l- Pi )l(l-p) A = r,/2d,.(l-A) Without loss of generality, let us start the cyclic schedule shown in Figure 3.1 at either Points PI or P5, and let us denote by i the index of the product we start with. According to Figure 3.1, the inventory of Product i behaves as shown in Figure 3.2. S t and s, are the maximum inventory and backlog levels respectively attained within a cycle of Product i. Note that s, is negative and S, is nonnegative. Let Q i = S, - s< . This quantity is known as the replenishment quantity in the inventory theory literature. Obviously, Q i must be positive. 33 Notice here that the production cycle is completely characterized when the quantities S ; , Q. and z i are determined. This, we do in the following. Based on Figure 3.2, it is not difficult to see that the total inventory holding cost and the total backlog cost over a cycle are given by 1 + S 2 Inventory cost = —c i 2 4(1 - A ) 1 s Backlog cost = — c, 2 ' 3(1 - A ) 2 ' 4(1 -a) The average setup, inventory holding and backlog cost per unit time of Product i is then given by FAS t ,Q t =— + — l ■ + — ~— — :£lL -- T 2Td i {l-p i ) 2T4(1"A) Si Inventory A % i >v i / \J T; + ti T - x t - U 1 Figure 3.2: Inventory behavior of Part Type i. 34 The total average setup, inventory holding and backlog cost per unit time is given by 2 f F(S l ,S 2 ,Q 1 ,Q 2 ) = 2 T 2 4(1- a) u:sf + cn^-Qp (3.1) In the following, we show how T and Q i are related to r, , the time spent producing Product i at the demand rate. First, it is not difficult to see that % = Qi Pi cUl-A)' r _ __L_ yC . fftzfll* (m-l)£?4 (m-1) tT(l-A) ' (1-p) Here, m = 2 (two products) . Substituting T in the expression of Q i , gives (3.2) From the demand satisfaction constraint Td { = r,d, + t { p ( , we get Q,.=d;(l- A )(r-z;.). (3.3) But T = V fa + r, + <5|). Substituting for t, and rearranging terms, we get (3.4) a = q, r- 1 + Zd-P^T " f 1 -^ o o J j-t (3.5) 4(1- a) Where q t = — £. q. can be seen as the replenishment quantity when the (1 -p) r/s are all zero. It is clear that the independent variables of the model are the S/s and the r ; 's. Letting S and r be the vectors which components are the S/s and the r/s respectively, the total average cost per unit time can be rewritten as follows: 35 l-mf 1 1 *^»*J- i j^q. — • (3-6) (m-l)tf4 (m-1) The minimization problem can be stated as follows: ( P ) Minimize F(S, r) Subject to Q = <?, i + Z(i -/>,)*}/£ - (1-pKAh *-i. •••> w; S,. £0, r,£0, Q,£0, i = l, ..., m. 3.2.2 Optimal solution Our solution approach consists of two steps. First, we obtain the optimal S/s as a function of the r/s through the quantities Q { . Second, we ignore the nonnegativity constraints and solve for the r/s implicitly through the quantities £>. If the obtained solution is feasible (i. e., satisfies the nonnegativity constraints), then it is optimal. If the nonnegativity constraints are violated, then in the optimal solution of the constrained problem, either one or both nonnegativity constraints will be effective. Hence, we must distinguish the three possible constrained solutions. That is, (r x = 0, r 2 > 0), (r, > 0, r 2 = 0), or ( x l = 0, r 2 = 0) . A detailed procedure for obtaining the feasible optimal solution will be provided later in the section. Setting the partial derivatives of F(S,t) with respect to S. to zero and rearranging terms gives S; = -r-^Q„ t = 1, ...» m. (3.7) Cr + C. 36 The second partial derivative of F(S, r) with respect to S, is always positive. Therefore, F(S, r) is strictly convex in S f and S* given above is a global minimum for F(S, r). Substituting S* in F(S, r) gives F(t)= k + aq: + aqi (3 . 8) QJd l+ QJd 2 -S Before we proceed with the solution for the r/s, let us prove the following proposition: Proposition 3.1 F(t) is strictly convex in ^ and r 2 . Proof. The Hessian matrix of F(t) is given by *=4 H n H 12 "21 ""22 Where H,, = Kb 2 +Aa 2 <(l-A) 2 r 2 +^ 2 ((1-A)^, +^) 2 5 It is clear that H u is positive. The determinant of His given by 4KAa 2 d 2 [a(l- A )r 1 -hb((l-p 2 )r 2 +^)] 2 det H = mo 4KA 2 b 2 c%[a{(l-p 1 )T 1 +^) + b(l-p 2 )r 2 ] 2 T 6 4<%<%A l A 2 a 2 b 2 [(l-p 1 ){l-p 2 )T 1 T 2 -{(l-p l )T l +S){(l-p 2 )T 2 +6)] 2 37 which is clearly positive. Hence, all the minors of H are positive, and H is positive definite. The result follows immediately. ■ As a consequence of Proposition 3.1, F(t) has a unique global minimum r*. In the following, we provide the optimal solution for the cases where (ri > 0, x 2 > 0), (zi = 0, z 2 > 0), (*i > 0, r 2 = 0) and the case (r x = 0, r 2 = 0). Cgse(r 1 >0,T 2 >0): Taking the partial derivatives with respect to r, and r 2 , setting to zero and rearranging terms gives the following system of two nonlinear equations: [A(2 A -l)Qf +^(l-2p 2 )<? 2 2 +2(p 2 (d 1 /d ! )A +(l-p 2 )&/d 1 )A 2 )Q 1 Q 2 -2*AAp& +A s d l {l-p 2 )Q 2 )-K = A 2 (2p l -l)Q 2 2 +A 1 (l-2p 1 )Qf + 2{p 1 (^/d l )A 2 +(l-p 1 ){cl 1 /^)A 1 )Q 1 Q 2 -2%A s &p l Q 2 +Adi{l-Pi)Qi)-K = Lemma 3.1: The linear equation A4Qi = A4»Qj ^ s a solution of the above system. proof: substituting Q 1 = (A 2 /i4 1 )(d 2 /d 1 )Q 2 in the first equation of (3.9) gives the second equation. ■ Solving (3.9) using Lemma 3.1, gives the following optimal values of Q l and Q 2 for the unconstrained problem (i.e., without the nonnegativity constraints) : or- ^Wi+^^r, (1-p) (i-A) (3.10a) 38 / <hAn Q u 2 = i + Ji +2 f (l-A) , (1-A) ri4 r 2 4> (3.10b) The superscript u stands for unconstrained. Note that, at this stage, we have no guarantee that i^ and i£ satisfy the nonnegativity constraints. Case ( r x = 0, r 2 unconstrained) : A similar procedure as in the previous case leads to the following expression for C^ and <? 2 \ ^-^(<(i- A f A +^) Case (r x unconstrained, r 2 =0): In a similar fashion, we obtain (if+AWf) (3.11a) (3.11b) (3.12a) (3.12b) Case (Tj = 0,r 2 =0) This case is trivial. Indeed, setting t x and r 2 to zero in (3.5) gives Q l =q l and<? 2 - q 2 . 39 Lemma 3.2: The following procedure gives the optimal solution of the constrained problem. Procedure 1 Compute g" and C% using (3.10a) and (3.10b); Compute t\ and t\ as follows: J..QE. °lh S; (3.13a) d\ (1-aK £«Sl ^ 8. (3.13b) d, (1-AWi JF % > and £ > THE7V zj = x\, t 2 = x\, Q\=Ql andQ' 2 =Q»; STOP ELSE Compute Q" and Q\ using (3.1 la) and (3.1 lb); Compute $ and t£ Using (3.13a) and (3.13b). IF t* > THEiV Calculate C 2 = F(<, z£) using (3.8); ELSE C 2 = oo. END Compute C\ u and C£ using (3.12a) and (3.12b); Compute r" and r 2 Using (3.13a) and (3.13b). IF x\ > THEJV Calculate C 1 = F«, z£) using (3.8); ELSE C, =oo. 40 END Let Q? = g, and £> 2 U = q 2 => «J = and t u 2 = 0; Calculate C = F(r", r 2 ) using (3.8); ((*ii *£M#.Q»)) = argmin{C , Cj.cJ; Proo/: expressing r, and r 2 as a function of Q" and Q 2 , gives T _Ql QTa ,. r, = d 4. (1-aK ra _ or cgft j, d, (l-^Jd, Now, if rj and r 2 are both nonnegative, then the solution of the unconstrained optimization problem is feasible and since F(t) is convex in T x and r 3 , it follows that the solution is optimal for the constrained problem as well. If the nonnegativity constraints are violated, then the optimal solution is obtained by comparing the costs of the three possible cases (effective constraints) and picking the minimum. For the cases where only one of the nonnegativity constraints is effective, if we obtain a nonfeasible solution (i.e., the other constraint must be effective too), we set its cost to infinity. This is because when only one nonnegativity constraint is effective, we set the variable for which the constraint is effective to zero and we solve an unconstrained optimization problem with respect to the other variable. Hence, we might obtain a nonfeasible solution. ■ 41 Once we obtain the optimal rj, r 2 , Q[, and <?*, we go back and calculate the optimal S f , s t , f, (i = 1,2) and T, which completely characterizes the optimal cyclic schedule. Points PJ(J=1,2, . . ., 8) of Figure 3.1 are given as follows: ( > Pl = P2 = P3 = P4 = P5 = S 2 - #A +a 2 p l d i /d 1 [l-p i ) S 2 - §& + Sjd^/d, (1 - a) - lid, fa-«w > V °2 y 'Si - S 3 d, + sjLjyJdt (1 - p 2 )] (3.14a) (3.14b) (3.14c) (3.14d) (3.14e) (3.141) (3.14g) (3.14h) Where s,, S,., Q, and r, (i » 1,2) are the optimal values calculated using Procedure 1, (3.7) and s, = S t - Q i . v j ( S, - 4d + s 2 djjjd 2 (1 - a ) - ^d ' P6 = °i - ^"i ■" V P7 = < S 2 J P8 = , S 2 - ^idz. 3.2.3 Numerical examples In this sub-section, we compare our results with other models found in the literature. We show that by allowing a certain production time at the demand rate, we achieve significant savings. 42 The examples used for the comparison are those proposed by Boctor (1982) to illustrate his algorithm for a two-product Basic Period schedule model. These two examples were then used by Lee and Surya (1989) to compare their new algorithm to Boctor's. Also, we use these two examples with the model of controllable production rates (no production at the demand rate) of Buzacott and Ozkarahan (1983). Notice that no backlog is allowed for the two examples. The case where backlog is not allowed is obtained by letting c. (the backlog cost rate) go to infinity. Which is basically equivalent to replacing y i by c* for i = 1,2 in the above formulas. In this case, the a^'s become all zero. Example 1 i d, u t 4 h < 1 20000 2 27000 160000 162000 0.0125 0.0250 15 25 0.005 0.004 Example 2 i 4 v t S; h «r 1 3500 2 46500 100000 • 100000 0.5 0.3 2500 18500 0.150 0.005 The optimal solution for the examples is as follows: Example 1 : x\= 0.4282, t\ = 0.4815, t[ = 0.0900, t' 2 = 0.1 1 1 1, T* = 1.1483. Example 2: rj = 6.0248, r 2 = 0.0000, t\ = 0.2521, t\ = 6.1509, T* = 13.2278. 43 Table 3.1 shows the computational results of this comparison. The upper part of each entry in Table 3.1 shows the average cost for each schedule, while the lower part shows how much savings are made with our schedule. For instance, in Example 2, our schedule average cost is 51.3% lower than that of Boctor's. Although we expected the Basic Period schedule to perform better than the Common Cycle schedule, the results show that by introducing a production at the demand rate, we achieve savings up to 66% on the average cost of the schedule. Notice that, even with a controllable production rate but without a production rate at the demand rate (Buzacott and Ozkarahan), we obtain lower average costs. Tab Example 1 Example 2 e 3.1: Comparison savings for the two-product model Our model 72.0 3403.8 Boctor 118.7 64.9% 5149.3 51.3% Lee and Surya 119.2 65.6% 4153.0 22.0% Buzacott and Ozkarahan 94.3 31.0% 4144.1 21.7% 3.2.4 Special cases As mentioned earlier, based on the optimal solution, we distinguish four different cases. These are shown in Figure 3.3 ( a, b, c and d). In this section, we derive some results that will be useful for Chapters 4 and 5. First, we derive a condition for which we have the case shown in Figure 3.3a. As will be shown, this corresponds to a system where setup costs are negligible and the system is heavily loaded. The transient solution of this case will be studied in Chapter 4. This case constitute the back bone of the optimal transient solution of the more general case studied in Chapter 5. Another 44 special case is obtained when the setup times are negligible but the setup costs are significant. We will show that in this case t\ and r 2 cannot be simultaneously zero. Figure 3.3a: r, = and r 2 = 0. Figure 3.3c: t x > and r 2 = 0. Figure 3.3b: r, = and r 2 > 0. Figure 3.3d: r, > and r 2 > 0. Figure 3.3: Possible shapes of the cyclic schedule 45 Case of zero setup costs In the following, we formally derive the condition for which the cyclic schedule in Figure 3.3a is optimal. This condition is given in Theorem 3.1. Theorem 3.1: Without loss of generality, let Part Type 2 be the part type such that Yi&i ^ TA- ^ th e setup costs fc, and fc, are zero, then r* = and r* = if and only p+p 1 -1^0. If p+p 1 - 1 < , then z\ = and r* = if and only if r2 d 2 / ri d 1 ^(l-A)/(l-/>-A)- Proof: The proof is based on the Karush-Kuhn-Tucker (KKT) optimality conditions. Substituting Q l and Q 2 in the objective function of Problem (P) above gives the following optimization problem. r t and r 2 are the decision variables: [jC + -fl r t (r +{ai -IK + a 2 vj + \h 2 {T +a 1 r 1 + {a 2 -l)r 2 f) Minimize F(t) = ± - - <- T +a 1 r, + a 2 r 2 Subject to t { >0, i = l,2 where H { = ^(l- A ),i«l,2; a,=(l- A )/(l- /3 ),i = l,2; T = S/(l-p). Since F(t) is strictly convex in r, and r 2 and the constraint domain is convex, it follows that the KKT optimality conditions are necessary and sufficient (which establishes the if and only if part of the theorem). Let g^i, ^ 2 ) = ~ T i (i = l>2). The KKT optimality conditions are given as follows: 46 ju f >0 for i = l,2 where r* is the optimal value of r, (i = l,2). VF(t 1} t 2 ), V^r^rJ and Vg 2 {T v T 2 ) are the gradients of F(r v t 2 ), &(«•„ r 2 ) and fir 2 (r 1} r 2 ) respectively. Letting r* = [i = 1,2) and calculating the gradients and substituting in the KKT conditions, gives: a-alH^^-^-^-O; 2 J * 2 r 2 " 2 ^ + f^_iV 2 -^-/, 2 = o. 2 U J ' TZ For the KKT conditions to hold true, we must have ft, > (i = 1,2). Hence, we must have r 2 1 2 J 2 r 2 2 1 2 J Substituting a,, a 2 , H x and H 2 by their expressions gives r 2 ((i- A)(/ 2 ^-M)+prM) . iv S — > 2 „ ^((p+ a-i^^+Ci-a)^) *~ 2 ' Letting K = , we get the following conditions: (1-pJrA-KO+prt^zO; ( a ) (p+A-i)rA+(i-A)rA^o. (&) 47 Notice that condition (a) is always satisfied since y 2 d^ > y^. For condition (b), if p + p 1 -l>0, then it is always satisfied. If p + p 1 -1<0, then it is satisfied only if y&ly&<{\- Pl )l{\-p-p x ).m The following corollary follows immediately. Corollary 3.1: Assuming that Part Type 2 is the part type such that y 2 d^ > y^ . If the setup costs hq and fc, are nonzero, then r* = and r* = if and only if o < k < (2?/2)niin{(i-AXyA-yA)+/^A(p+A-i)yA+(i-A)rA}- The following theorem shows that, when the setup costs are zero for both part types, at most one part type will be produced at the demand rate, during the cyclic schedule. Theorem 3.2: If the setup costs fc, and k^ are zero, then at least one of the r/s will be equal to zero. Proof: We prove the theorem by contradiction. Assume that t\ > and r* >0, then Q[ and Q\ are given by (3.10a) and (3.10b), where K = 0. In this case, we have <?; = <? = 2qAr» 2<&dift Using (3.13a) and (3.13b), we get . ZckrA 2q 1 y 2 d 2 p 1 J>Q . d 2 {a 2 y 1 d l + a^d, ) d, (1 - a )( a^ + a^d, ) 48 . 2q 1 y 2 d 1 2q 2 y 1 d^ 2 ^ n r 2 = -. r -. r - O > U . &\<wA + wAi) AQ- PiK^rA + <hrA) Substituting q x and q 2 by their expressions, simplifying and rearranging terms gives: (i-A)rA-Q-+Pi)rA >0 ^d 0--f\)rA-Q>+P2)rA>o. ft4 (1+a) M (i-a) Which cannot be. Hence, in the case of no setup costs, r\ and t 2 cannot be simultaneously nonzero. ■ Case of zero setup times The steady state solution for this case can be obtained using the same procedure used for the case where setup times are nonzero. The following theorem shows that in the case of zero setup times, the cyclic schedule must have at least one segment with production at the demand rate (i.e., it is not possible to have the case of Figure 3.3a). Theorem 3.3: If the setup time is zero for both part types, then there is at least one segment within the cyclic schedule corresponding to a production at the demand rate. That is, if ^ = <J 2 = 0, then it is not possible to have i = r 2 =o. Proof: Assuming that t\ = and r* = 0, and letting £-> in the objective function introduced in the proof of Theorem 3.1, gives: ( | HX | H 2 T^ F(0,0) = lim- _K_ T n '0 Since T = 5/(1- p) and K>0. 49 Hence, this solution cannot be optimal. In the case of zero setup times, the general structure of the cyclic schedule is as shown in Figure 3.4. The coordinates of the cyclic schedule in this case are given as follows: Pl = s a +s l (d I M) A /(i- A ) (3.15a) P2 = P3 = P4 = JS 2 +s,(dsM)A/(l-A)- tA \ S 2j (3.15b) (3.15c) Figure 3.4: General structure of the cyclic schedule in the case of no setup times. 50 P5 = P6 = f (3.15d) ^ + s 2 (d 1 /^)p 2 /{l-p 2 )- rA] (3.15e) J PI = P8 = \ S 2j (3.151) where s it S,, Q, and r, (i = 1,2) are the optimal values calculated using the same procedure as in the case of nonzero setup times. 3.3 The Multi-Part-Tvpe Case In this section, we solve the multi-part-type problem numerically. We adopt Zoutendijk's Algorithm (see Bazaraa and Shetty, 1979). To simplify the algorithm, we exploit the structure of the problem and implement the algorithm using closed form expressions. First, we give the mathematical formulation of the multi-part -type problem. 3.3.1 Problem formulation It is not difficult to extend the two-part-type problem formulation to the m-part-type case. Using the expression of T given by (3.4) and substituting Q i (i = 1,2) given by (3.5) in the objective function given by (3.8) gives the following formulation: P ) Minimize F(t) = if JT +j^H t {T - T t f ) Subject to r, > 0, i = l, ..., m 51 where H { = rA[l-p { ), t»l, -., m; «, =(l-A)/( 1 -/°)» * = 1 > •••> m ' T = d/(l-p). T is the length of the cyclic schedule when the r,'s are all zero. Notice that the above formulation is expressed as a function of the t t 's only. It is not difficult to show that S,, given by (3.7), is still optimal, which completely characterizes the optimal cyclic schedule. Before we proceed with the algorithm, we propose the following conjecture. Conjecture 3.1: Based on the strict convexity of the two-part-type objective function, we conjecture that the m-part-type objective function is strictly convex in the r/s. The following proposition guarantees that the optimal r,'s cannot be infinite. Proposition 3.2: Problem (P) has always a finite optimal solution. proof: Assume that the optimal solution of Problem (P) is not finite and label the part types so that r t , r 2 , ..., r n are infinite and r n+1 , r n+2 , ..., r m are finite (n<m). Without loss of generality, let X = r, = r 2 = ... = r n . Substituting and rearranging terms in the expression of T gives It is clear that if X -> », then t lf r 2 r n will all go to infinity. Now, substituting in the objective function, we get 52 F(r) = f 1 1 ( \ ^ + i2 2 y n H I (M 1 +M 2 -l) 2 +-^ 2 y m H, M.+M^^l where M. = — + T ""a, and M, = V , " m a, — . limM 2 = 0. Hence, limF(r) = oo, which is obviously not optimal since we can always choose a finite feasible solution with a finite objective value. ■ Since the objective function is strictly convex and the solution is always finite, it follows that Problem (P) has a unique global solution. 3.3.2 Algorithm First, we state Zoutendijk's algorithm (Bazaraa and Shetty (1979)) for rninimizing a differentiable function / in the presence of linear constraints of the form Ax < b . Initialization step: start with a feasible solution x 1 . Let te=l, and go to main step. Main step: 1. given x k , suppose that A 7 and b T are decomposed into (Al,Al) and (b T ,b T ) so that A 1 x k = b x and A,x fc < b 2 . Let d k be an optimal solution to the following problem. (PI) Minimize V/(xJ r d Subject to >ijd<0; -l<dj <1, j = l, ..., m, if Vf{x k fd k = 0, stop; x k is the optimal solution. Otherwise, go to Step 2. 2. let X k be an optimal solution to the following line search problem: 53 (P2) Minimize f(x k +kd k ) Subject to 0<k<k B where k =< max minlfa./d, : d, > 0} if not alld, are negative oo if alld, are negative b = b 2 -A 2 x k ; d = A,d k . Let x k+1 = x k + k k d k , identify the new set of tight constraints at x k+1 and update \ and A, accordingly. Let fc:= fc + 1, and repeat step 1. Notice that Problem (PI) is a linear programming problem. Since in our case A = -I m and b = 0, where 7 m is the mxm identity matrix. (PI) can be solved by inspection as follows. Let Vi^(r) be the I th component of the gradient of the objective function F(r) in problem (P). v^)>^' r - r ' l ; H ' (r - f,) -* r( ' l) . 13.10 Let 3 k = {i : if = 0} be the set of indices for which r, is equal to zero at iteration k (tight constraints). Then (PI) can be rewritten as follows: Minimize j£ ?%[?)&, Subject to < d,. < 1 , i e 3 fc ; -l<d, <1, t"«*3 fc . Since for a linear programming problem the optimal solution is always at a vertex, the following procedure gives the optimal solution: Procedure 2 IFie3 k THEN 54 IFVF i (f)>0 THEN d,:=0; ELSE d ; :=l; END ELSE IFVF i (f)>0 THEN d,:=-l; ELSE d,:=l; END END The optimization problem (P2) is a line search problem which can be eliminated by calculating the optimal value of X analytically. Differentiating F(t k + Xd k ) with respect to k and setting the result to zero gives the following quadratic equation in X . ^ETH.M-d,)^ + 2M 2 £;i 1 m H,.(M 1 -d i ) 2 A + 2M 2 £ i TH 1 .(M 1 -d 1 .)(M 2 - I f) = o, where M, = j£ a,d, and M 2 =T + j£ arf . The positive solution of the above quadratic is given by 2K+Y m K M | Z^-i i M M 1 M i I ZZ H iM-*f (3.17) 55 b and d can be obtained as follows. Since b = and A, = -I n (loose constraints), it follows that b { = rf and d. = -df for i'f«3 r The algorithm can now be stated as follows: Initialization step: start with a feasible solution r 1 . Let k: = 1 , and go to main step. Main step: 1 . Identify the set 3 k ; Apply Procedure 3.2: -> d k ; IF Vf(ffd k =0 THEN STOP; r* is the optimal solution. ELSE GOTO Step 2. 2. £>,.:= if and d,:= d* for ig3 fc . imin^./d,. : d, > 0} if not alld,. are negative oo if alld,. are negative £:= nun (A*,A m J; T /c+1 :=^+l fc d fc ; k:=k + l; GOTO Step 1. 3.3.3 Numerical example with backlog allowed As an example, consider the ten-product problem proposed by Bomberger (1966), which we extend to the case where backlog is allowed. The backlog cost is 30 times that of inventory holding cost. For this example, we 56 have normalized the data as follows: the demand rate of each part type is set to 1, the maximum production rate for each part type is divided by its corresponding demand rate, the inventory cost rate of each part type is multiplied by the corresponding demand rate and finally, the backlog cost rate of each part type is multiplied by the corresponding demand rate. The data for the example is as follows. d; u t 1 1 15.3 0.500 130 0.20896 6.2688 2 L 23.5 0.750 200 0.03188 0.9564 3 100.0 0.500 110 0.02321 0.6963 4 L 18.8 0.125 10 0.01667 0.5001 5 : L 47.5 0.250 30 0.01063 0.3189 6 ] L 80.0 0.125 20 0.00490 0.1470 7 ] L 400.0 1.000 310 0.00375 0.1125 8 ] L 300.0 0.250 50 0.00223 0.0669 9 ] L 150.0 0.125 5 0.00170 0.0510 10 ] 300.0 0.125 5 0.00027 0.0081 The optimal solution is as follows: i i t; s; * 1 109.78 1.74 24.1 -0.80 2 5.80 126.4 -4.21 3 1.36 130.7 -4.36 4 7.28 125.0 -4.17 5 2.87 129.3 -4.31 6 1.71 130.4 -4.35 7 0.34 131.7 -4.39 8 0.45 131.6 -4.39 9 0.91 131.2 -4.37 10 0.45 131.6 -4.39 57 T* = 136.46 days Average cost per day = $13.05. 3.4 Summary In this chapter, we have studied the steady state version of the Production and Setup Scheduling Problem. This turned out to be a special Economic Lot Scheduling Problem with controllable production rates and backlog is permitted. We derived the optimal solution of the two-product problem in closed form. For the multi-product problem, we proposed a very simplified version of Zoutendijk's algorithm which requires neither solving a Linear Programming sub-problem for finding a feasible direction, nor performing a line search procedure to determine the next improving solution. Comparison with previous results reported in the literature revealed that savings up to 66% can be obtained when controllable production rates are allowed. CHAPTER 4 OPTIMAL TRANSIENT SOLUTION OF A TWO-PART-TYPE MANUFACTURING SYSTEM WITH SETUP TIMES AND NO SETUP COSTS DETERMINISTIC SYSTEM 4.1 Introduction In this chapter, we study the transient behavior of the deterministic two- part-type system. We consider the case where setup costs are negligible and the system is heavily loaded. In Section 4.2, we present the optimal location of the cyclic schedule in x-space. In Section 4.3, we show how we partition the x- space into two mutually exclusive major regions 9t" and 91°. In Section 4.4, we determine the optimal transient solution in Region 9f . In Section 4.5, we determine the optimal transient solution in Region SR U , using an algorithmic approach. In Section 4.6, we analytically determine the optimal transient solution in Region 9T. We conclude this chapter with a summary of results. 4.2 Optimal Location of the Cyclic Schedule We assume that the setup costs fc, and fc, are zero for both part types and that the conditions of Theorem 3.1 of Chapter 3 are fully satisfied. This guarantees that the optimal cyclic schedule has the shape shown in Figure 4.1. That is, there is no production at the demand rate at the axis. In this case, The location of the cyclic schedule is obtained by minimizing the average inventory and backlog costs incurred during one complete cycle and is given by the following points in x-space: 58 59 d/ x 2 c \ ^^> A \ / X1 B Figure 4.1: Optimal cyclic schedule A = U ,B = (B [B 2 1 Lc- yC* ,D = fa K D 2j ,i = fi \ \^2j (4.1) where _cf_(l-A) (cT+cDfl-p)^ « + q) (1-/?) < (1-A) (c 2 + + c-)(l-p)^' B 2 = - ^- -^^; (c 2 + + c-)(l-p) c c 2 - (1-A) ^ (c 2 + +c-) (1-a^ A=" « + cT)(1-p) <Sc^, A = (1-A) (c 2 + + c 2 )(l-p) &£> - ^ti, ; (4.2) (4.3) (4.4) (4.5) 60 Here and elsewhere in the dissertation, p i = d i /U i (i = 1,2) denotes the utilization factor of the machine by Type i parts, and p = p t +p 2 is the utilization factor of the machine. For the problem to be feasible, we need to have p<l. S= S 1 +S 2 is the sum of the setup times. The total average cost over one cyclic schedule is given by Average Cost = 1 g£ Ezfl) ^ + 1 c ^ UzBl g*. 2(cT + c 1 + )(l- /7 )^ 2(c-+c 2 + )(l- / ?p 4.3 Partition of the x-Space into Two Mutually Exclusive Regions Based on Definitions 2.3 and 2.4 (Chapter 2), the optimal transient solution corresponds to a state trajectory emanating in x-space and reaching the cyclic schedule with minimum total cost and in finite time. The problem consists in finding such optimal trajectories for every initial production surplus level in x-space. Before we start deriving the optimal transient solution, we notice the following facts. Fact 4.1: Since the cost coefficients in the HJB equation of Theorem 2.1 are function of x, given an optimal choice er*(r), t< r<t s , the jc-space can be partitioned into mutually exclusive regions (Gershwin (1993)), each corresponding to an optimal setup state c? = (o\,cf 2 ,o\ 2 ,<J 2l ) and an optimal production rate u* = (u^ ,iC,) which belongs to the set D.*, defined in Chapter 2. Another important fact drawn from the HJB equation is given as follows: Fact 4.2: given an optimal choice cr*(r), t < T<t s , J T , the transient cost component does not depend explicitly on time and is piecewise quadratic. Therefore, the cost coefficients of the HJB equation in Theorem 2.1 are linear in 61 x. This implies that the boundaries of the regions in x-space must be linear (Gershwin, 1993). For the subsequent analysis, the following important quantities need to be defined. First we need to find analytic expressions for the equations of line [A,D) (LI), line (B,Q (L2), line parallel to (A,D) containing Point C (L21) and the line parallel to {B,C\ containing Point A (LI 2) (Figure 4.2). These line equations are easy to calculate since we know the coordinates of the points A, B, C, and D (from the steady state optimal solution) . Line Ll:d^(* I -4) + 4(l- A )(x 2 -A ! )«0; (4.6) Line L2:d 2 (l-p 2 )(x l -C i ) + d 1 p 2 (x 2 -C 2 ) = 0. (4.7) LineL12: d^(l- p 2 ){ Xl - A i ) + d 1 p 2 (x 2 - A 2 ) = 0; (4.8) LineL21:c3 2 p 1 (x 1 -C 1 ) + d,(l- A )(A: 2 -C 2 ) = 0; (4.9) Having the equations of Lines LI and L2, we can easily find the coordinates of the intersection point I, which is also of importance in the subsequent analysis. Intersection Point 7 is given by Ii" (1-A) (I-P) cf+c," ■d l 8-d i S 2 h = (1-A) 0--P) V C 2 +C 2 ■dj-dA (4.10) Notice that intersection Point I can be located in any quadrant of the production surplus space, depending on the values of the penalty function coefficients and the setup time duration. Without loss of generality we index the parts such that Part Type 1 is the part type with the larger setup time (i.e., S 1 > 8 2 ). To obtain the optimal transient solution, we divide the surplus space into two mutually exclusive major regions: 5R U and 9f . 62 L21 LI ^s. D " *2 \ \ A \ A x, V 2 V L12 Figure 4.2: Illustration of lines LI, L2, L12 and L21. 9t u is the region in x-space located below Lines L12 and L21. 9t° is the region in x-space located above Lines LI 2 and L21 (Figure 4.3). Algebraically, * u -«*,*») I d i (\-p 2 )(x 1 -A l ) + d l p 2 (x 2 -A 2 )<0; d^ix, -q)+^(i-AK«8-c a ) < o}. <R° = {(x 1 ,x 2 )|d 2 (l-p 2 )K-A) + d 1 P 2 (^-A 2 )>0 c4a(«i -CJ + ^I-aKx, -C 2 ) > 0}. We divide the analysis into two parts. In the first part, we determine the optimal transient solution for all initial surplus levels in Region 91°. In the second part, we do the same for all initial surplus levels in Region 91" . 63 L21 / *2 iVAa \ /\ *i M u v \ NL12 Figure 4.3: Regions M" and 91° 4.4 Optimal Transient Solution in Region 9*° As mentioned above, our task is to bring the initial surplus levels in x- space to a point on the cyclic schedule. Because, once on the cyclic schedule, the production surplus will behave according to it thereafter (since the machine is reliable), and we know that operating according to the cyclic schedule is optimal (from the steady state solution). The solution approach for initial surplus levels in Region 5R° is obtained by inspection. However, it is easy to see that the solution is indeed optimal. For instance, from Theorem 2.1 and Fact 4.1, we know that the optimal production rates are piecewise constant and are selected from a very small set. Hence, from the system dynamics (differential equations (2.1)), it is easy to verify that the optimal trajectories in x-space are 64 piecewise linear and therefore it would be not difficult to see that they are indeed optimal. We will provide a detailed solution for the case where Intersection Point I is located in the first quadrant (Figure 4.3). The positive location of Intersection Point / may be the most common in practice. It has been reported in the literature that the backlog and inventory costs, cT and c* (i = 1,2), are usually chosen at least in the ratio of 5 to 1 (because backlog usually results in sales and customer's good will losses and other undesirable effects). The setup times S t (i = 1,2) are of the same order of magnitude. Hence, from the expression of the coordinates of Intersection Point /in (4.10), we have: </(c; + c+) > Sj/^ + S 2 ) for i = 1,2; j = 1,2 and j * i . Hence, i, and i 2 are positive numbers. In any case, the solution procedure is exactly the same for Intersection Point J located in any other quadrant in x-space. Throughout the dissertation, we assume that initially the machine is not setup to either part type. For initial surplus levels in Region W, we can state our problem as follows: Subject to constraints (2.1)-(2.6), find an optimal trajectory (i.e., with minimum cost) emanating from an initial point x(t) eW in x- space and reaching the cyclic schedule at a point x(tj, where it would be possible to move according to it thereafter. The solution to the above problem depends on the initial surplus levels. Therefore, we further partition Region 91° into mutually exclusive regions G, Gl, and G2 by introducing two linear boundaries LG1 and LG2. Line LG1 (LG2) is the line for which the surplus level of Part Type 1 (Part Type 2) is exactly equal to d l S 1 (d^). This partition is shown in Figure 4.4. Its justification will be apparent from the solution below. 65 4.4.1 Initial surplus levels in Region G In this case, both surplus levels are positive. Therefore, the production is ahead of the demand for both parts and it is optimal to let the surplus levels deplete at the maximum possible rate. That is, we set the production rates to zero. However, we also want to reach the cyclic schedule as quickly as possible. Hence, we need to decide which setup we must prepare the machine for, and when is the best time to do it. To answer these questions, we further partition Region G into six regions: Gil, G12, G21, G22, HI and H2 located in the first quadrant (Figure 4.5). Gl L21 LG1 Figure 4.4: Further partition of Region 9T. 67 containing Point /with direction (d,,d,). The optimal trajectory is obtained by setting the production rates to zero. This corresponds to a trajectory moving downward in the southwest direction with speed (-d,,-^). When the trajectory touches the boundary of Region G12 on Line LI 2, we immediately start a setup for Part Type 2. At the end of the setup change the trajectory touches the cyclic schedule at a point on the Segment [I,C\ (trajectory r 12 in Figure 4.5). Initial surplus in Region G2 1 Region G21 is defined as follows: G21 = {(x,,x 2 )\-d^(x i -A l ) + d i (x 2 -A 2 )>0; d i (x t -I i )-d 1 {x a -I 2 )>Q; d 2 p 1 (x 1 -C,) + d 1 (l- A )(x 2 -C 2 )^0}, where the first constraint corresponds to the half space above Line LA. The second constraint corresponds to the half space below Line LI. The third constraint corresponds to the half space above Line L21. LA is the line containing Point A with direction (dpdj. In this case, the optimal trajectory is obtained by first setting the production rates to zero. This corresponds to a trajectory moving downward in the southwest direction with speed (-d,,-^). When the trajectory touches the boundary of Region G21 on Line L21, we immediately start a setup for Part Type 1 . At the end of the setup, the trajectory touches the cyclic schedule at a point on the Segment [I,A] (trajectory r 21 in Figure 4.5). Initial surplus in Region G22 Region G22 is defined as follows: G22 = {( Xl ,x 2 )\ x 2 - S& > 0; djx, -A,) -<L,(x 2 -A,) > 0}, 68 LG2 A Figure 4.5: Optimal trajectories emanating at points in Regions Gil, G12, G21, G22, HI and H2. where the first constraint corresponds to the half space above Line LG2. The second constraint correspond to the half space below Line LA. In this case, the optimal trajectory is obtained by first setting the production rates to zero. This corresponds to a trajectory moving downward in the southwest direction with speed (-dj,-^). When the surplus of Part Type 2 reaches level S 2 d, (that is the trajectory hits Line LG2), we immediately start a setup for Part Type 2. At the end of the setup, the surplus level of Part Type 2 is exactly zero, that of Part Type 1 is still positive and we still have not reached the cyclic schedule. Therefore, we need to decrease the surplus level of Part Type 1 and keep the surplus level of Part Type 2 at the zero level until we touch the cyclic schedule. To do this, we produce Part Type 2 at the demand rate d, . This corresponds to 69 a trajectory moving along the x, axis with speed (-d^O), until the cyclic schedule is touched at Point B' (trajectory r 22 in Figure 4.5). Initial surplus levels in Region HI Region HI is defined by m = {(x 1 ,x 2 )\ x l -8 l d l >0; d 2 (x l -C l )-d i [x 2 -C 2 )>Q; d 2 p,(x,-C l ) + d l (\-p 1 ){x 2 -C 2 )>0; -d l (\-p 2 )(x l -A l )-d i p 2 (x 2 -A 2 )>0}. The first constraint of this set corresponds to the half space below Line LC. The second constraint corresponds to the half space to the right of Line LG1. The third constraint corresponds to the half space above Line L21 and the fourth constraint corresponds to the half space below Line LI 2. Region HI can be further partitioned into two mutually exclusive regions: Hll and H12. Region HI 1 is defined by Hll = {(*,, x 2 )| J^-^d, >0 d 2 (x 1 -C 1 )-d 1 (x 2 -C 2 )>0; -<*,(*, - 9 n ) + d 1 (x 2 - g 12 ) Z 0; -d 2 {l-p 2 ){x,-A x )-d l p 2 {x 2 -A 2 )>0}, where g 1 = (g n ,g 12 ) is the point in x-space given by g n = S& and g l2 = C 2 + (C, - ^djd^/djl - Pl ) . Point g 1 is the intersection of Lines L21 and Lgl (Figure 4.6). For initial surplus levels in Region Hll, the optimal trajectories are obtained in the same manner as for initial surplus levels in Region Gl 1 (trajectory r u in Figure 4.6). Region H12 is defined by H12 = {(x,,x 2 )|d 2 p 1 (x 1 -C 1 )-d 1 (l- A )(x 2 -C 2 )>0; 70 -^ta - 0ii) +d i(*2 - sy ^ °; For initial surplus levels in Region H12, the optimal trajectories are obtained in the same manner as for initial surplus levels in Region G21 (trajectory r 12 in Figure 4.6). Initial surplus levels in Region H2 Region H2 is defined by H2 = {(x 1 ,x 2 )\x 2 -S 2 d i >0; d 2 (x 1 -I 1 )-d 1 {x 2 -I 2 )>0; ^(1 " Pa)(Xi ~ A) " 4A(*a " A) * 0; -d 2 p 1 (x 1 -C 1 )-d 1 {l- Pl )(x 2 -C 2 )>0; -d^-AH^te-AJX)}. The first constraint of this set corresponds to the half space above Line L12, the second constraint corresponds to the half space above Line LG2, the third constraint corresponds to the half space below Line LI, the fourth constraint corresponds to the half space below Line L21, and the fifth constraint corresponds to the half space above Line LA. Region H2 can be further partitioned into two mutually exclusive regions: H21 and H22. Region H2 1 is defined by H21 = {{x i ,x 2 )\d 2 (x 1 -IJ-d l (x 2 -I 2 )>0; -d 2 p 1 (x 1 -C 1 )-d 1 (l- Pl )(x 2 -C 2 )>0; -d 2 (x 1 -g 21 ) + d y (x 2 -g 22 )>0), where g 2 = (g 21 ,g 22 ) is the point in x-space given by 02i = A + (A ~ <W«W48(1 ~Pi) and &2 = <W 71 Point g 2 is the intersection of Lines LI 2 and Lg2 (Figure 4.6). For initial surplus levels in Region H21, the optimal trajectories are obtained in the same manner as for initial surplus levels in Region G12 (trajectory r 12 in Figure 4.6). Region H22 is defined by H22 = {(x lt x 2 )\x 2 - S& £ 0; -d i (x 1 -A 1 )-d 1 {x 2 -A i )>0; -d i p 1 (x 1 -C 1 )-d 1 (l-p 1 )(x 2 -C 2 )>0; -d>(x 1 -92i) + d i( x 2-922)>°}- For initial surplus levels in Region H22, the optimal trajectories are obtained in the same manner as for initial surplus levels in Region G22 (trajectory r 22 in Figure 4.6). Figure 4.6: Partition of Regions HI and H2. 72 4.4.2 Initial surplus levels in Region Gl Region Gl is defined by Gl = {(x 1 ,x 2 )\-x l +S 1 d l > 0; d 2 p 1 (x 1 -C 1 ) + d 1 (l~p 1 )(x 2 -C 2 )>0}. The first constraint of this set corresponds to the half space to the left of Line LG1 (Figure 4.7). The second constraint corresponds to the half space above Line L2 1 . Notice that in this region, the surplus level of Part Type 2 is always positive and the surplus level of Part Type 1 is less than ^d, . Therefore, we need to produce Part Type 1 so that we can reach the cyclic schedule. Suppose that the surplus level of Part Type 1 is positive. In this case even if we started a setup for Part Type 1, we would endup with a backlog of Part Type 1, since during the setup for Part Type 1, 8^ amount of Part Type 1 is depleted. The optimal trajectories emanating from Region Gl and leading to the cyclic schedule are obtained as follows: First we start a setup for Part Type 1. This corresponds to a trajectory moving downward with speed (-d,,-^). When the setup is completed, x^ is negative and x 2 is positive. Hence, we produce Part Type 1 at the maximum machine production rate to eliminate the backlog of Part Type 1 . This corresponds to a trajectory moving southeast with speed (U t -(2,,-dj). When the surplus level of Part Type 1 becomes zero, that of Part Type 2 is still positive (Point where the trajectory hits the x 2 axis), and we still have not reached the cyclic schedule. Therefore, we need to decrease the surplus level of Part Type 2 and keep the surplus level of Part Type 1 at the zero level until we touch the cyclic schedule. The optimal way to do this is to produce Part Type 1 at the demand rate dy. This corresponds to a trajectory moving downward along the x 2 axis with speed (0,-dj) until the cyclic schedule is touched at Point D' (trajectories r, and rj in Figure 4.7). 73 4.4.3 Initial surplus levels in Region G2 Region G2 is defined by G2 = {(x v x 2 )\-x 2 + S& > 0;d,(l - a)K " A) + <& 2 (*2 - A) * 0}. The first constraint of this set corresponds to the half space above Line LG2 (Figure 4.7). The second constraint corresponds to the half space above Line L12. In this region, the surplus level of Part Type 1 is always positive and the surplus level of Part Type 2 is less than 5^. This is the symmetric case of initial surplus levels in Region Gl. That is, the optimal way to get to the cyclic schedule is to setup the machine for Part Type 2 first, produce Part Type 2 at the maximum machine production rate until its surplus level becomes zero. At that point change the level of production of Part Type 2 to d> and continue producing this part until the trajectory touches the cyclic schedule at Point B'. (trajectories r 2 and if 2 in Figure 4.7 ). To summarize the control actions in Region 91°, let x = (a,b) be the vector of initial surplus levels. Then, we have • If x eGl: 1 : Setup the machine for Part Typ e U 2: Produce Part Type 1 at the rate U 1 until the surplus level of Part Type 1 becomes 0; 3: Change the production rate to d,; 4: When the surplus level of Part Type 2 becomes x 2 = A, +A 1 d 2 pJd 1 (l -f\), switch to the control actions of the cyclic schedule. • If xeGllUHll: 1 : Do not produce either part type; 74 LG2 *~^ G2 Figure 4.7: Optimal trajectories emanating from Regions Gl and G2. 2: When the surplus level of Part Type 1 becomes £&, start a setup for Part Type 1 ; 3: Produce Part Type 1 at the demand rate d,; 4: When the surplus level of Part Type 2 becomes x 2 = Aj +.A,d 2 /3,/d l (l -#), switch to the control actions of the cyclic schedule. • If jceG12UH21: 1 : Do not produce either part type; 2: When the surplus level of Part Type 2 reaches Level l 2 , immediately start a setup for Part Type 2. l 2 is given by i 2 =(l- A )(b-(d 2 /d 1 )a) + (d 2 /d 1 )(l-p 2 )A+P 2 A; 75 3: When the setup is over, switch to the cyclic schedule control actions. •If xeG2lUH12: 1 : Do not produce either part type; 2: When the surplus level of Part Type 1 reaches Level l x , immediately start a setup for Part Type 1 . ^ is given by l 1 =p 1 (b-(d 2 /d l )a) + (d 2 /d 1 )p 1 C 1 + (l-p 1 )C 2 ; 3: when the setup is over, switch to the cyclic schedule control actions. • If xeG22UH22: 1 : Do not produce either part type; 2: When the surplus level of Part Type 2 becomes 6^, start a setup for Part Type 2; 3: when the setup is over, produce Part Type 2 at the demand rate d^; 4: When the surplus level of Part Type 1 becomes x 1 = C i +C 2 d 1 /7 2 /d 2 (l -/7 2 ), switch to the control actions of the cyclic schedule. • If xeG2: 1 : Setup the machine for Part Type 2; 2: when the setup is over, produce Part Type 2 at the rate U 2 ; 3: When the surplus level of Part Type 2 becomes 0, change the production rate to c^ ; 4: When the surplus level of Part Type 1 becomes x, = Cj + C 2 d 1 p 2 /d 2 (l- p 2 ), switch to the control actions of the cyclic schedule. 76 Notice that all optimal trajectories, emanating in Region 5R°, reach the cyclic schedule by moving downward and that only one setup for either part type is performed before the cyclic schedule is reached. In other words, for initial surplus levels in Region 9T, the cyclic schedule is always reached from above with only one setup. 4.5 Optimal Transient Solution in Region 9T (Algorithmic Approach) In this section, we develop an algorithm to obtain the optimal trajectories emanating in Region 9t u and reaching the cyclic schedule in finite time. Before we proceed with the algorithm, we establish the following facts and definitions. In the previous subsection, we showed that for initial surplus levels in Region 9T, the cyclic schedule is always reached in finite time from above with only one setup. For initial surplus levels in Region 9V , if we immediately started a setup for a part type, at the end of the setup, we would miss either Segment [A,D] if we set up for Part Type 1, or Segment [B,C\ if we set up for Part Type 2 (which are the target segments for the surplus to stay on the cyclic schedule) . To bring the surplus levels to either segments in finite time, we need to generate a surplus excess for the part type the machine is set up for, so that when we switch to the production of the other part type, we endup on or above the appropriate segment of the cyclic schedule. After the setup, we need to produce the part type we set up the machine for. Based on Theorem 2.1, we can set the production rate to zero, to the demand rate, or to the maximum rate. It is clear that if we set the production rate to zero, both surplus levels deplete and the generated trajectory moves in the direction (-d,,-d,) and hence further drifts away from the cyclic schedule (Figure 4.8). If we produced at the 77 demand rate the part type the machine is set up for, we would keep the surplus level of this part type constant, while the surplus of the other part type depletes. In this case, the generated trajectory will move parallel to one of the axis in x-space in the direction that further drifts away from the cyclic schedule (Figure 4.8). If we produced at the maximum production rate the part type the machine is set up for, the surplus level of this part type increases, that of the other part type decreases. This way it is possible to generate an excess surplus for the part type being produced. The following fact formally states this result. Fact 4.3: For initial surplus levels in Region 9T , the only way to progress toward the cyclic schedule is by producing at the maximum production rate whenever it is possible. Based on Fact 4.3, the trajectories emanating in Region 5R" either move along direction (-d,,-^) during a setup, along direction {U t -d^-dj) during the production of Part Type 1 or along direction (-d i ,U 2 -d^) during the production of Part Type 2. Therefore, the cyclic schedule cannot be reached from below (i.e., from points in Region 9T) in finite time. To be able to reach the cyclic schedule in finite time, starting in Region 9t", we should bring the surplus levels in Region 9t° , then apply the optimal control actions of that Region (since we know that those controls lead to the cyclic schedule in finite time) . Doing so will generate a surplus excess of either part type. To reach the cyclic schedule with the least amount of excess, we must bring the surplus levels to the boundary of Region 9T (Line Li/' i = 1,2; j = 1,2 ; i * j) with minimum cost. Once on the boundary Li/', we switch to Part Type j and produce this part type (this corresponds to a trajectory moving along Line L/) until the cyclic schedule is reached at Point D if j=l or Point B if j=2 (Figure 4.8). The following fact formally states this result. 78 *! / \ (-4,^r4)\ D"^ (-4,0) —er—-y IvsA *i / / / / / (-4,-4) k / / / \ /\ \ / \ / s \ Y (^i-4,-4)W yL12 (-4,-4) i \ i \ i \ (0,-4) \ \ \ \ N /\ NL2 ^ Figure 4.8: Illustration of Facts 4.3 and 4.4. Fact 4.4: To reach the cyclic schedule in finite time, starting with initial surplus levels in Region 9t", the trajectory leading to the cyclic schedule must touch the boundary Lij of Region 91° just before switching to Part Type j and reaching the cyclic schedule at one of the points BorDon the cyclic schedule. Definition 4.1: We say that a trajectory is following Direction D t , if it moves parallel to Line Li (Figure 4.2) in the direction of increasing x { . That is, the machine is producing Part Type t (i = 1,2). Remark: Since for surplus levels in Region SR" the machine always produces at its maximum rate, the trajectories will move either along Direction Dj , along Direction D 2 or along Direction D (during a setup for a part type) . 79 Definition 4.2: We call a D^n-step trajectory (i = l,2; n>l), a trajectory that performs alternately m i setups-production runs of Part Type i and rrij setups-production runs of Part Type j (/>i), with the initial setup for Part Type i and the last segment touching the cyclic schedule at Point B or Point D. If n is even then, m, = m } = n/2. If n is odd then, m i = (n + l)/2 and rrij =(n- 1)/2. Figure 4.9 shows a D 1 -2- step and a D 2 - 2 - step trajectories. We can now state the problem as follows: Given surplus levels in Region 9t u , find a D,. -n*-step trajectory that minimizes the cost of reaching the cyclic schedule and satisfying constraints (2.1)-(2.6). Where, n* is the optimal number of steps (setups-production runs) before reaching the cyclic schedule and D t , is the first direction to follow (i.e., the initial setup). /\^L21 LI / \^S. 1^,2-stepV D "\ ^\ \ ^ 2 ^A^A L 12 \L2 \ Figure 4.9. D 1 -2- step and D 2 - 2 - step trajectories. 80 As a consequence of Facts 4.3 and 4.4, the optimal trajectory leading to the cyclic schedule in finite time is at least D ( -2- step . Since, we need a setup change at the initial point and one more, once the boundary of Region 9T is reached. This leads to the following proposition. Proposition 4.1: Let C"(x)(i = 1,2) be the cost of moving along a JD, - n - step trajectory with initial surplus point x in surplus space. An upper bound on the optimal cost of reaching the cyclic schedule is given by min^Cfix)}. Proof: Assume that the optimal trajectory with initial point x is D ( . - n * -step (n*>2). We know that an optimal trajectory originating at Point x in surplus space is at least D ( - 2 - step . It follows that, C£» < C?(x) and Cg{x) < C 2 2 [x). Therefore, C?(x) < rmn{C*(x)}. ■ Based on Facts 4.3 and 4.4, The optimal trajectory emanating at a point in Region 9T and leading to the cyclic schedule in finite time can be obtained as follows: Given an initial surplus point in Region 9T, we choose the first setup and calculate the cost of the trajectory leading to the cyclic schedule with two setups only. At this point, we have a D t -2- step trajectory. Where, i is the initial setup for Part Type i. The next step is to try to lower the cost of the current trajectory by introducing a setup for the other part type so that the cyclic schedule is reached at the opposite side. If the cost can be reduced, the obtained new trajectory is a D, - 3 - step trajectory. We keep trying to reduce the cost of the current trajectory by introducing, each time, a setup before the cyclic schedule is reached until we cannot lower the cost anymore. The resulting trajectory is an optimal D { -n-step trajectory (provided we start with a setup for Part Type i] emanating in Region 5R U and reaching the cyclic 81 schedule in finite time. In the same manner, we obtain the optimal D j -n- step 0«] trajectory starting with a setup for Part type j first. The optimal trajectory would be the one with the lower cost. The numerical solution of various examples suggests that the above procedure be further simplified based on the following conjecture. Conjecture 4.1: The initial setup of the optimal trajectory is given by i* = argmin 1 . w {C I 2 (x)}. Conjecture 4.1 implies that we only need to compare two D t -2-step trajectories instead of comparing two D t -n~ step trajectories. Which is simpler and faster to compute. In other words, the initial setup is optimal and once we obtain the initial setup, we do not need to go back and check whether starting with the other setup is optimal. Although we do not have a mathematical proof for this observation, the numerical solution of an extensive number of cases agree with Conjecture 4.1. Conjecture 4.1 along with the above procedure lead to the following algorithm which is called Direction Sweeping Algorithm (DSA) . Direction Sweeping Algorithm Notation : x(0): Initial surplus in x-space (Input); (x(0),x(l), ..., x(k)}: Optimal trajectory (Output); J x : Optimal trajectory cost (Output); C(x(k)): Cost of the trajectory up to the point x(k); Cgm(Y,Z): Cost along Segment [Y,ZJ. Dr: Direction of search. Dr = D 1 or Dr = D 2 as defined in Definition 4.3. Dr = 1 ifDr = 2 and Dr =2 ifDr = 1 . 82 C£, (X) : Cost ofDr-n-step trajectory starting at the point X in x-space. I(X,Dr): Intersection point of the line containing X and with direction Dr, with Line L^—, where L^.— = L 12 or L^-^ = L 21 , depending on the direction of search Dr. P^.-P^B andP 2 =D. Main Algorithm : // Initialization: m:=l; STOP:=FALSE; // Find an upper bound on the cost of the optimal trajectory emanating from x(0): (J x ,Dry.= min{C%(x(0))y, // Update trajectory: k:=k+l; x[k):= xik-^ + S^-cL,-^); Dr:= not Dr; // Update trajectory cost up to the point x(k): C(x(k)):= Cgm{x(k - l),x(fc)); WHILE STOP * TRUE DO swap(Dr,Dr); // Add a new step to the trajectory: {Cl(Y(m)),Y(m)):= min {Cl{Y) + Cgm(x(k),Y,Dr~)}; IF C(x(k)) + Cgm(x(k), Y(m)) + C% (Y(m)) < J x THEN // Update trajectory: k:=k+l; 83 x(k):=Y(m); k:=k+l; xfkJr-xflc-lJ + ^H^-d,)! // Update upper bound on the optimal cost: J x := C(x(k)) + Cgm(x(k),Y(m)) + Cl{Y(m)); // Update trajectory cost up to the point x(k): C(x{k)):= C(x(k)) + Cgm(x(k - 2),x(k - 1)) +Cgm(x(k - l),x(k)); ELSE // Cannot add another step to the current trajectory. // Then the trajectory x(k) is optimal. STOP .— TRUE; END END DO swap(Dr,Dr); k:=k+l; x(k):=I(x(k),Dr); k:=k+l; xikY-xik-V + S-i-^-d,); k:=k+l; x(k):=P^; //Optimal trajectory cost: C(x(k)):=J x ; To help understand the DSA algorithm, we illustrate it step by step using the following example. Throughout the example, the reader is referred to Figure 4.10. Given an initial surplus vector x{0), we go through the following steps: 84 Step 1: Compute C?{x{0)) and C£(jc(0)). Say that C 2 (x(0)) < C?{x(0)), and let J x be the upper bound on the optimal trajectory cost and the initial direction of search is Direction 2. That is, J x = C 2 (x{0j) and Dr=2. Now, we know that the optimal trajectory is at least D 2 -2-step. Then, we try to add another step to the trajectory and see if we can lower the upper bound J x . This we do in the next step. Step 2: for all Y e[x(l),/(x(l),2)], find the best D l -2- step trajectory that starts from Y and leads to the cyclic schedule at Point B. To do this, we perform a line search along Segment [x(l),/(x(l),2)], using the golden section algorithm (see Bazaraa and Shetty (1979)). The optimal D x -2 -step trajectory starting at an initial point on Segment [x(l),4x(l),2)] is obtained for Point Y[l) with cost Cf(Y(l)). Now, we check whether it is optimal to add this D 1 -2-step trajectory to the trajectory up to Point x(l) and obtain a D 2 -3-step trajectory with a lower cost. The cost of this new D 2 - 3 - step trajectory is Cl(x(0)) = C(x[l)) + Cgm(x{l),Y(l)) + Cf(Y(l)). Where C(x(l)) is the cost of the trajectory up to Point x(l) starting at Point x(0), and Cgm{x{l),Y[l)) is the cost of Segment jx(l). V|])l. Now, suppose that C 2 (x(0))<J x . Then, it is better to move along this new D 3 - 3 - step trajectory. The new upper bound on the optimal trajectory cost is J x = C 2 (x(0}}. In the next step we check whether it is possible to add another step to the current trajectory and reach the cyclic schedule at Point D. Step 3: Let x(2)=Y[l) and x(3) = x(2) + S 1 (-d 1 ,-d i ). Now, the situation is identically the same as if we started with Point x(l). Therefore, we need to go through Step 2 once again. Suppose, we did that and found a new upper bound J x = Cl (x(0)). Then, it is an improvement to add a new step to the trajectory and obtain a new trajectory D 2 - 4 - step . 85 IM1),2) % \I(x(3),l) ^^XlfxfO),!, Figure 4.10: Illustration of the DSA algorithm. Step 4: Let x(4)=Y[2) and x[5) = x(3) + <X 2 (-d 1 ,-d 2 ). Similarly, we go through Step 2. Suppose, we did that and could not lower the upper bound on the optimal trajectory cost. Hence D 2 -4-step is optimal. To complete the trajectory, we let x(6)=l(x{5),2), x(7) = x{6) + S 1 [-d 1 ,-d i ], and x(8)=D. The optimal trajectory is given by {x(0),x(l ),..., x(8)} and its cost is given by J x . Notice that each time we try to add a new step to the current trajectory, we sweep all possible trajectories in the direction of search and pick the one that minimizes the cost of reaching the cyclic schedule (this is done by means of the line search procedure). This is why the algorithm is called Direction Sweeping Algorithm. In the following, we prove the validity of the line search and the uniqueness of its solution. For this, we use the algorithm notation. 86 Proposition 4.2: The cost function C 2 ^ (x) is strictly convex in x. Where, x e[x(k),I(x(k),Dr)]. Proof: Let Cl(x) = Cgm(x,x' ,0) + Cgm(x' ,I(x' ,Dr),Dr) + Cgm(I(x' ,Dr),x" ,0) + Cgm(x n ,P^,Dr). Where, Cgm(X,Y,l\ is the cost of Segment \X,Y\ when direction I is followed. Here, direction is the South-West direction. In other words, the machine is being set up. In general, we can write the cost along Segment [X,Y\ in x-space as follows: Cgm(X,Y,l) = f j- ** {yf - x? } where X = (x,,x 2 ), Y = (y lf y 2 ) and h, = c* , if x,. and y i are positive and h t = -c~, if x. and y, are negative for i=l,2. If 7is fixed, it is easy to show Cgm\X,Y,1) is a strictly convex function of X. In particular, Cgm(x",P^,Dr) is a strictly convex function of x" for fixed P^ (recall that P^ is either Point B or Point D on the cyclic schedule). Since x" can be expressed as a linear function of P^, Cgm(I(x' ,Dr),x" ,0) + Cgm(x" ,P-,Dr) is a strictly convex function of ^x",Dr) for fixed P- . Continuing the reasoning in the same manner, and since the domain [x(fc),7(x(fc),£>r)] is convex, it follows that C^(x) is a strictly convex function of x e[x(k),I(x(k),Dr)]. ■ Using Proposition 4.2, it follows that the line search in the algorithm always provides a unique optimal solution. Theorem 4.2: The DSA algorithm gives a unique optimal trajectory with the optimal cost C!f (x(0)), where x(0) is the initial point of the trajectory in x- space, m* is the optimal number of steps in the trajectory and i* is the initial direction to follow. 87 Proof: To see that the trajectory obtained by the DSA is indeed optimal, consider an initial point x(0). Then, the DSA algorithm can be described by the following recursive expressions: C*(x(0)) = min{q 2 (x(0))}; 1-1,2 C^ +1 (x(0)) = min{c^(x(0)); min {C(Y) + Cl{Y)}} y^jc{2m-3|J(j;(2m-3)j„)] where m=2,3,...,m*; i, = i *, i 2 = i*, i 3 = i*, ...; i* is the starting direction of search, C™(x(0)) is the cost of the Dj*-m-step trajectory and C(Y] is the cost of the current trajectory up to point Y. The remaining quantities are the same as defined by the algorithm. These recurrence relation and Proposition 4.2 guarantee that each step added to the trajectory is optimal. Therefore, by the principle of optimality, the trajectory obtained is optimal. Furthermore, by Proposition 4.2 the optimal trajectory is unique. ■ Figure 4.11 shows an optimal trajectory starting at Point (-10,-10) for the following problem data: d, = 2, d> = 3, U^ = 6, U 2 = 6, S 1 = 2, 8 2 = l, cf = 10, c 2 = 10, c," = 50, and cj =50. Notice that, the optimal trajectory is D 1 - 5 - step. Figure 4.12 shows different optimal trajectories starting at different initial points. Notice that, the setup switching policy is a kind of corridor policy. The walls of the corridors have negative slopes with respect to the axis. Also, the corridor is linear. Which means that the boundaries in Region 9T are linear as stated in Fact 4.2. 88 -40 -20 20 Production Surplus of Part Type 1 40 Figure 4.1 1: Optimal trajectory obtained by the DSA algorithm. -25 -15 -10 -5 Production Surplus of Part Type 1 10 Figure 4.12: Linearity of the boundaries in Region 9T. 89 4.6 Optimal Transient Solution in Region 9i" (Analytical Approach) In this section, we derive the optimal transient solution in Region 91" analytically, based on the DSA algorithm of the previous section. Together with the optimal solution of Region 9T, we obtain the complete analytical optimal transient solution. The algorithmic approach can be used when one tries to extend the problem to the multi-part-type case. The analytical approach is always desirable especially in the real-time scheduling context. It usually involves the use of a formula, which is an excellent tool for studying extreme cases where one is interested in the behavior of the solution with respect to one or many parameters of the problem. In addition, it provides insights to the more general case and indicates the pertinent parameters of the system. 4.6.1 Preliminary discussion and results According to Fact 4.4, to reach the cyclic schedule in finite time, we need to perform a setup on one of the two boundaries LI 2 or L21. However, starting from a point in SR", we may need more than one setup to bring the surplus trajectory to the cyclic schedule. For instance, as illustrated in Figure 4.13, starting at Point in 9t" with an initial setup change for Part Type 2, we have a choice to change the setup at Point 1 on the boundary L21 and reach the cyclic schedule at Point D, or have setup at Point 2 and then at Point 3 on L12 and reach the cyclic schedule at Point B. If the cost can be reduced, we may want to initiate a setup on the segment [2',3] and swing back to boundary L21. As a consequence, the optimal trajectory leading to the cyclic schedule in finite time will be at least D i -2- step, since we need a setup at the initial point and one more, once the boundary of Region 91° is reached. 90 "\1 2 L \ c \ A *i \ ° 0' B\ \ \ \L12 3 \ \ Figure 4.13: Illustration of the setup switching policy Figure 4.12 of the previous section shows the application of the DSA algorithm to an example where optimal trajectories starting at different initial points are shown. Notice that the setup switching policy is a special corridor policy (see Sharifnia et al. (1991)) with two windows. Based on numerical results, we observe that the general optimal setup switching policy is as illustrated in Figure 4.14. In this case, the surplus trajectory will bounce back and forth between the two walls of the corridor {<M and <B2 in Figure 4.14) until the trajectory passes through one of the windows of the corridor. When this happens, the trajectory will encounter either Line L12 or Line L21 then bounces back and follows either Line L2 or Line LI and reaches the cyclic schedule at Point D or Point B (Figure 4.14). Therefore, the optimal setup switching policy is completely characterized when the corridor walls <B1 and <B2 91 and the corridor windows defined by Points Wl, W2, VI and V2 in x-space are determined. The purpose of this section is to obtain, analytically, the optimal setup switching policy in Region 9T . This is achieved by determining the equations of the boundaries <B1 and <B2 as well as Points Wl, W2, VI and V2. One of the advantages of the analytical solution is that, once we have the equations of Boundaries (Bl and <B2, we eliminate the line search procedure in the DSA algorithm and substitute it with an intersection point calculation, which is clearly much faster and much simpler to implement (as a formula). Before we procedure any further, we state some preliminary results that are needed in the subsequent analyses. L21 *2 <B1 ~--~-^^yi L\"~~—-. yC \I Corridor Window A / 7\ / x, \ Corridor Wall Surplus Trajectory \ W Av2 1 \ \ * / A \ \ ^ 1 \ \ N / \ \ W \ \ y ^ \ A \ \ / x v \ / v M / N M / s \ / N \ / \ V / v |\ / v l v * \ A N V \ x \ 1 v \ \^~4 \ \L12 \ 1/ \ \ x \ )^ \ v v \ f 2 \<B2 Figure 4.14: Illustration of the setup switching policy. 92 We denote by Line L(Y), the line going through a point Y = (y 17 y 2 ) in x- space and parallel to Line L2. We denote by Line L'(Y), the line going through a point Y = (y 1 ,y 2 ) in x-space and parallel to Line LI. The equations of L(Y) and L'(Y) are given by LineL(Y): ^(1- p 2 )(x 1 -y 1 ) + d i p 2 (x 2 -y 2 ) = 0; (4.11) LineL'(Y): ^(x, - yj + c^l- A )(x 2 -y 2 ) = 0. (4.12) Let Z 1 = (z 11 ,z 12 ) be the intersection point of L(Y) with L21 and Z 2 = (z 21 ,z 22 ) be the intersection point of L'(Y) and LI 2. Then, Z x and Z 2 are given as follows: (i-p) I (1-a) J * 12 = ^^ j-(*M)«fc -tt^t% +(*M)flfi +(1-A)C 2 1; (4.13b) (l-/o) I (1-P 2 ) J - = i^|--=^»i-(^MW a +(i-A)A+(^/4i)AA}; ( 4 - 14a ) (l-p) l (l-A) J z,, = z 22 - TT^i^M^ +d- Afefe -(d,M)AA -^rl- ( 4 - 14b > We denote by C XY the cost of Segment [X,Y\ in x-space, when X and V are within the same quadrant in x-space. For X = (x lf x 2 ) and Y = (y lf y 2 ), we have C xy = f fil fof -xf} (4.15) where h, = c* if x, and y, are positive and h, = -ef if x, and y. are negative for {-1,2. Remark: If X and 7 are not in the same quadrant, we can always break the segment into two or three parts each contained within one quadrant in x- 93 space. Then, the total [X, Y] segment cost is the sum of the two or three segment costs. 4.6.2 Equations of Boundaries <B1 and <B2 In this subsection, we determine the locus of the points where it might be optimal to initiate a setup for the other part type. In other words, how far should we go with the production of the current part type before we switch to the other. Notice that, at this point, we don't know whether this switch will lower the cost of the current trajectory or not. This will be decided later on, when we determine Points Wl, W2, VI and V2. Indeed, if the switch happens to be within the corridor windows, then it is optimal not to switch and continue producing the current part type until either Boundary L12 or Boundary L21 is encountered, where it is optimal to switch to the other part type and reach the cyclic schedule at Point B or Point D. To determine the equations representing the corridor walls in jc-space, we start with an initial surplus point in the third quadrant of the x-space. It will be shown that the equations of the corridor walls (Boundaries <B1 and <S2) do not depend on the initial surplus. To be able to follow the derivation, the reader is referred to Figure 4.15. First, we derive the equation of Boundary ®1. Before we start the derivation, we need to determine the coordinates of Points X = (x 1 ,x 2 ) and V = (v 1 ,v 2 ), given X ={x° + S 2 d l ,x° +S&) (without loss of generality). In this case, X belongs to Line h{X' Q ) and V is the intersection of Line L12 and Line h'(X), where X' = (x?,x°). Then, using (4.11), (4.14a) and (4.14b), we have: 94 Figure 4.15: Corridor walls x 2 =(d 2 /d 1 )^-^i(x 1 -x 1 ) + x 2 ; A (l-/>) I (1-a) J (1-p) [ U-AJJ As mentioned before, the objective is to determine the location in x-space of the switching point X that minimizes the cost of the trajectory emanating at Point X in x-space and reaching the cyclic schedule at Point B (Figure 4.15). Notice that the current trajectory in Figure 4.15 is a D 2 -2-step trajectory given by {X^X'^X^X'^DJ} . Now, if we introduce a setup switch at Point X, the new trajectory will be a D 2 - 3 - step trajectory given by {X ,X' ,X,X',V, V',B,I) . 95 At this stage, we have no guarantee that the new trajectory will have a lower cost than the current one. This will be decided later when we determine Wl, W2, VI and V2. All we are saying at this point, is that if switching to the other part type would lower the cost of reaching the cyclic schedule, then the optimal switch must happen at Point X in x-space. Using the notation of Proposition 4.1, C 2 (X ) would be the cost of Trajectory {X ,X' ,X,X' ,V,V ,B,I}. Using (4.15), we have 2^(1 -a) 2cyi-/? 2 ) __L [tf + SAf q- +I?c*}+ (X2 ° )2 ° 2 — (x, +SAf<$ +—-£ -ltd, 2d/^ 2 ^'^ ^ J 2d 2 (l-p 2 ) 2d^* ™ 2d,(l- A ) 2"2" Notice that, if we substitute x 2 , l?„ t> 2 by their expressions as a function of Xj, C 2 (X ) will be a quadratic function of x,. Hence, C 2 (X ) may have a unique minimum. A necessary condition for Xj to be a mimmum of C*(X ) (and sufficient in this case) is given by dC 3 3 (X ) Q Differentiating x 2 , i^, v 2 with respect to x,, we get dXj p 2 Now, differentiating C 2 (X ) with respect to Xj and expressing v, and v 2 as a function of Xj and x 2 and setting the result to 0, gives the following equation of a line in x-space: 96 (Bl: M& + M 2 x 2 +M 3 = (4. 16) where m, = -^—c,-- flg; < - wmwa)^ 1 ^* ( 4 - 17a ) (i-a) (i-a)(i-p) l 1 -/ 7 ) (1-/7) p 2 P 2 (l-/?) A^-^q + MA). ( 4 - 17c ) Notice that <3i contains Point C on the cyclic schedule. This result is not a surprise since, if we started with an initial surplus X belonging to Line L12, Point C on the cyclic schedule would realize the minimum for C 2 (X ) and therefore is a candidate for initiating a setup. With a similar procedure, we obtain the equation of Boundary <B2, where it is possible to switch from Part Type 1 to Part Type 2. (B2: N& + N 2 x 2 + N 3 =0 (4.18) where y ,-[ <tM )-fl- < -lq-- ' t -ft'' 1 -ft' q- i (4.19a) (1-p) a A(W) ^-j—c-- *g - <-(4MWa)^q; (4.1%) (1-A,) (1-A-)(1-P) ( 1_ P) iV 3 = -(^A+^ 2 A)- ( 4 - 19c ) It is worth mentioning that in the previous section, we observed from numerical solutions that the optimal setup switching boundaries in x-space are linear. It is verified here by the equations of the corridor walls. 97 4.6.3 Corridor windows In this section, we determine the corridor windows, which are characterized by the coordinates of Points Wl = (w u ,w 12 ), W2=(w 21 ,u> 22 ), Vl=(v ll} v 12 ) and V2= (v 2l ,v 22 ) in x-space (Figure 4.14). Wl and W2 are the points on the corridor boundaries where an additional setup does not reduce the cost of reaching the cyclic schedule. Geometrically, this means that beyond these two points Wl and W2 on the boundaries Gil and <B2, it would not be optimal to initiate a setup (that is, the cost of the current trajectory cannot be lowered). First, we derive the coordinates of Wl and W2. The coordinates of VI and V2 are obtained as a function of those of Wl and W2. We derive the coordinates of Wl first. To be able to follow the derivation the reader is referred to Figure 4.16. As we mentioned above, Wl is the point in x-space that realizes the equilibrium between the cost of Trajectory {W1,V1,V1',) and the cost of Trajectory {W1,W1',Z1,Z1',I}. Wl must be a point of <B1 since we know that <B1 is the set of optimal locations where we switch to Part Type 2. Figure 4.16 shows an example of three trajectories t lt r 2 and r 3 . For trajectory r,, it is optimal to switch to Part Type 1. For trajectory r 2 , we are indifferent between switching to Part Type 1 or continuing the production of Part Type 2, switching setup at Point VI and then reaching the cyclic schedule at Point D. In this case the two trajectories have the same cost. For trajectory r 3 , it is not optimal to switch to Part Type 1 and it is optimal to continue the production of Part Type 2 and reach the cyclic schedule at Point D. Before we start the derivation, we need the coordinates of Point Z\ = (z u ,z l2 ) shown in Figure 4.16. In this case, VI belongs to L(W7)flL21 and Zl belongs to L'(Wl')f]L12. Once we have VI and Zl, VI' and ZV are easily obtained. VI and Zl are defined as follows: 98 Figure 4.16: Corridor windows 2 n = Wu + niJ ** = 7 ?i2 w i2 + n 2 ; y i2 = Xl2" ; i2+^l 2 where fti =7~t{WM)(M j /M 2 K1- A )- A }; (1-/?) (1-AK1-A) . . (1-Aift (W) A+- (1-P) -(d,/d 2 )(^ 2 +(M 3 /M 2 ))- (1-P) 4A; (4.20a) (4.20b) (4.21a) (4.21b) (4.22a) (4.22b) (4.23a) 99 H. " ^^{-A^MJA "7^7 A ~ Q.-p l )[MjM a )-(d 2 /d l )S 1 A (4.23b) r 12 = tt^tUW/MJ- (1 -pJteM)}; (l-yf) «u - -TT^Ci -^fWcy(C 2 +(M 3 /M 2 )); (I-/7) (l-p) /j,(l-^) , , (1-A)(1-P 2 ) ^ 1 AA , M /M * (l-p) (l-p) (i-p) (4.24a) (4.24b) (4.25a) (4.25b) Denote by TC U , the cost of Trajectory {W1,V1,V1',I} and by TC 12 the cost of Trajectory {W1,W1',Z1,Z1',I}. Using (4.15), we have TC, if tufa I& Jfci-WcT+Afoil | [<&-p 2 wl 2 c + 2 } d, " d, + " d,(l- A ) da(l-P 2 ) (4.26) 12 2 d, 4(1 "A) u£< [(gi 2 -^) 2 c 2 +^ 2 c 2 + ]' d, ^(1-pJ (4.27) As we mentioned above, Wl = (w n ,w 12 ) is the point in x-space for which we have TC U = TC l2 . Furthermore, Wl belongs to the boundary ®i. Therefore, Wl should satisfy the following two conditions: [rc il K 1 )-7U 12 (u; 11 ) = 0; \w 12 = -(MJM 2 )w 11 -(MJM 2 ). (4.28) 100 Notice that, substituting v llf v l2 , z n , z 12 and tu 12 by their expressions as a function of w n , makes TC U and TC l2 depend on w u only. After manipulation, we obtain the following system of equations: \w 12 = -(MjM 2 )w n -(M 3 /M 2 ) where (4.29) [A - iK" - 4i< , (r? 2 - i M JM 2 ) 2 )c; - fe t (430a) 2d,(l- A ) 2c^(l-p 2 ) 4(i-A) M 2 J 4,(1 -A) (4.30b) a _ MiiJMn ~ 2^)q - ( vfi - A 2 )q + 24(1 -A) 04 - (M 3 /M 2 ) 2 - 7 2 2 )c 2 + - ( v 12 - ^fci (4.30c) 2*4(1 -A) For most practical situations (finite c* and cr (z = 1,2), it can be shown that Wl is located in the second quadrant of the x-space. Hence, w u in (4.29) must be the negative root. In a similar fashion, we determine the coordinates of the Point W2 in x- space, by solving the following system of two equations: [#"4 +fi 2 w 22 +& = 0; where 101 n _ [lii ~ [NJNJ )< - n 2 j\ | [y 22 -l)c 2 - T] 22 c; , (432a) 2cMl- A ) 2d^(l- A ) ' A- cf- %i( v «-$<W d,(l-A) A 2^(1- A ) , Attl^ ~ 24^2 ) C 2 ~ ( ^2 ~ J 2 ) C 2 + . (4 32c ) 2^(1 -a) (I-P) %. = ^-{(d./d 1 )(N 2 /N l )(l-p 2 )-p 2 }; (4.33b) (1-p) v = ilz^L^^c i -A(^M)C 2 -(l-A)(^ 3 /iV 1 )-(d l /ci 2 )4d 2 }; (4.34a) (1-A I (1-a) J v 22 - ^-^(d.MKQ HNjN ^-f 1 -^ ^ + /^; (4.34b) r a , = 7 7^-(A(^/iV 1 )-(l-A)(d 1 /d 2 )}; (4.35a) (1-A r 22 - ^{(i-Al-A^MKW)}; ( 4 - 35b ) (i - p) (i-/?) (i-p) iy 22 in (4.31) must be the negative root. 102 4.6.4 The complete analytical solution At this stage, we have the complete analytical solution of the optimal production and setup control policy in Region 5R U . For initial surplus x = (a,b) in Region 5R U , the following procedure gives the optimal trajectory of the surplus levels in x-space. Procedure 1 STEP 0: *,(()):= a„X 2 (0):=b; fc-1; IF C?(X(0)) > C 2 2 (X(0)) THEN X^k^X^-S.cL; GOTO STEP 1; ELSE X 2 {k):=X 2 {0)-SA\ GOTO STEP 2; ENDIF; STEP 1: fc:=fc+l; XJfc):=(d^ 2 (M 3 +X 2 (fc-l)MJ + d 2 (l-/ ? JX^-l)Mj/(^(l-/ 72 )M 2 -d^ 2 M 1 ); X 2 (k): = -(MjMJXJk) - (MjM 2 ); !FX 1 {k)iW n THEN X 1 (fc):= ril X 1 (/c-l) + // 11 ; X 2 (fc):=r 12 X 2 (fc-l) + // 12 ; fc-fc+1; 103 x 1 {k^x 1 (k-i)-s 1 d 1 ; X a (fc):=X,(fc-l)-<5&; fc-fc+1; X(k):=I; GOTO STEP 3; ELSE GOTO STEP 2; ENDIF; STEP 2: k:=k+\; X 2 (fc):=(d> 1 (iV 3 +XJfc-l)iVJ + d 1 (l- A )X 1 (fc-l)^)/(d 1 (l- A )iV 1 -cV7 1 iV 2 ); X 1 [k):=-(NjN 1 )X 2 (k)-(NjN i ); IF X 2 {k) > w 22 THEN X l (k):=Y 21 X 1 (k-l) + M2l ; X 2 (k):=y 22 X 2 (k-l) + M22 ; fc:=fc+l; X^mXtlk-D-SA; k:=k+l; X(k):=I; GOTO STEP 3; ELSE GOTO STEP 1; ENDIF; STEP 3: OUTPUT. OPTIMAL TRAJECTORY = {X(0), X(l), ■■■■ , X(k-1), X(k)}; 104 Given the optimal trajectory of the surplus levels in x-space, it is not difficult to translate it into control actions. In the following, we give two examples which are solved using Procedure 1. Example 1 illustrates a nonsymmetric case where product 2 is more expensive than product 1. For this example, we start with a backlog of 10 units of product 1 and with a zero inventory level of product 2. Example 2 illustrates a symmetric case where both products have the same importance. For this case, we start with a zero inventory level for both products. In each case, the optimal trajectory is given as a collection of points in x-space (points where the trajectory changes direction). The points shown in italic are the points at which the trajectory follows the cyclic schedule thereafter. Example 1 Input: Cl + = 0.75, c; = 7.5, E7, ■ 7, 4 - 3, $ - 1; c 2 + = 1.25, c; = 18.75, U a m 5, d, = 2, 8 2 - 1.5. Initial surplus J5fl(0)=(-10,0). Output: Optimal Trajectory ={ (-10,0), (-14.5,-3.0), (-52.8,35.3), (-55.8,33.3) , (24.0,-6.6), (19.5,-9.6), (-13.9,23.7), (-17.8,21.7), (-2.3,14,4)=D,(22. 7, 1.9)=A,(18.2,-1.1)=B, (0. 7,16.4)=C, (-2.3,14.4)=D, ...}. Initial setup change to Part Type 2. Example 2 Input: q + = 1.0, q = 30.0, U l = 6,d 1 = 2,S l = 1; c 2 + = 1.0, c~ = 30.0, U 2 - 6, dj - 2, S 2 - 1. 105 Initial surplus X(0)=(0,0). Output: Optimal Trajectory ={ (0,0), (-2.0,-2.0), (16.3,-14.2), (14.3,-16.2), (-5.7,13.8), (-7.7,15.8), (12.6,-1.8), (10.6,-3.8), (-0.6,13.11), (-2.6,11.1), (-0A,9.6)=D, (11.6,1.6)=A, (9.6,-0A)=B, (1.6,U.6)=C, (-0.4,9.6)=D, ...}. Initial setup change to Part Type 1 . In the following section, we study some extreme cases which are of either practical or theoretical importance. Since we have the analytical solution, it is relatively easy to study such cases. 4.6.5 Special cases In this section, four special cases are investigated. Case 1 : No backlog is allowed for Part Type 2 (without loss of generality). Case 2: No backlog is allowed for both parts. Case 3: No inventory is allowed for part type 2. Case 4: No inventory is allowed for both part types. 4.6.5.1 No backlog for Part Type 2 This case can be encountered in practice where a backlog for a Part Type is undesirable either because of customer requirements or because of a very high penalty incurred when the products are not delivered in time. To obtain the optimal policy in this case, we simply let the backlog instantaneous cost of Part Type 2 go to infinity. That is c~ -» oo. As before, we determine the corridor wall first, then we determine the windows of the corridors (if they exist) . Corridor wall <B1: dividing the equation of <B1 by c~ in (4.16) and taking the limit when c~ -» <», we get 106 limM./c" = -KMWA)(l-p 2 )/(l-A limM 2 /c- = -{l-A)(l-/n,)/p a (l-P); limM 3 /c- = (d l /d 2 )(^ 1 / /?2 )((l- A )/(l-^))C 1 +((1- A )(l-p 2 )/ A (l-p))C 2 . Rearranging terms in the expression of <&1, gives the following line equation: d 7 p 1 (x 1 -C 1 ) + d 1 (l-p 1 )(x 2 -C 2 ) = 0. Which is exactly the equation of Line L21. Hence, when backlog is not allowed for Part Type 2, Boundary <B1 becomes Line L21. Corridor wall <B2: dividing the equation of (82 by c~ in (4.18) and taking the limit when c 2 -» oo, we get KmNjc; = 0; limJV 2 /c 2 =l/(l- A ); KkiNJc 2 = -6&l(l- p 2 ), which gives the following line equation for <B2: This result was intuitively expected. Indeed, x 2 cannot be negative. Hence, the switch to Part Type 2 must occur exactly when its surplus reaches Level S 2 (L, , so that after the setup is completed, the surplus level of Part Type 2 is not negative. Corridor Windows: To determine the corridor window, we divide expression (4.30a), (4.30b) and (4.30c) by c 2 then take the limit when c 2 -> oo and get: lim^/c- = a; = -((MJM 2 fc + 2 + jfi 2 )/2dv(l-p 2 ); lima 2 /c 2 = a; = -((MjM 2 )(MjM 2 )c + 2 + /&,( v 12 - <y a d,))/d,(l- P 2 Y> 107 1 \ \ y^\ \ \ h k jn/ X\ \ \ / D ®2 \\ \ 2 *2 yc K_4 V V V *L B Figure 4.17: Case where c~ -» oo. lima 3 /c- = < = -((M 3 /M 2 fc 2 + +( v 12 - S&f )/ 2^1 -p 2 ). Now, if we compute the discriminant in the quadratic o^w\ x + a 2 w u + a™ = 0, we get Kf - 4«)«) = -4c+ 2 [{MjM 2 )( v 12 - ^dj + (MJM,)^] 2 , which is obviously negative. This means that the cost of reaching the cyclic schedule from Point D, is always less than that of reaching the Emit cycle from Point B. Hence, the corridor window is the whole corridor wall given by the equation of Line L21. Figure 4.17 illustrates this policy. 4.6.5.2 No backlog for either part type In this case, we have c," -> oo and c~ -> oo. In a similar fashion as the previous case, we can show that, <B1: x, = §■& and <32: x 2 = S^. This means 108 that Region 5R U becomes "forbidden". That is, all initial surplus levels must be in Region 9T (otherwise the cost would be infinite) and therefore, the optimal control actions of Region 91° apply (see Figure 4.18). Figure 4.18: Case where cj" -» oo and c~ -» oo. 4.6.5.3 No inventory for Part Type 2 This case happens in practice when there is no space available for a product and the latter has to be delivered to the customer at once after being produced. Another situation is when the inventory cost for the product is much higher than its backlog cost. In this situation it is better to deliver the product late than store it for a while before being delivered to the customer. This is usually known as "produce to order system" in the literature. To obtain the 109 optimal policy in this case, we simply let the inventory instantaneous cost of Part Type 2 go to infinity. That is c 2 -»■ ». Corridor wall <B1: dividing the equation of <B1 by c 2 in (4.16) and taking the limit when c* -> », we get KmMjc; = 0; limM 2 /c 2 + = -l/p 2 ; limM 3 /c 2 + = 0, which gives x 2 = 0. This result was intuitively expected. Indeed, x 2 cannot be positive, hence the switch to Part Type 2 must occur when the backlog of that part type is completely eliminated. Corridor wall <B2: dividing the equation of <B2 by c 2 in (4.18) and taking the limit when c 2 -» oo, gives limiV 1 /c 2 + = -(4 2 M)A/(l-/'); c 2 - * 00 limiV 2 /c 2 + = -AP 2 /(1-A)(1-P); limiV 3 /c 2 + = (d 2 /d 1 )U/(l-/ 7 ))A+(^ 2 /(l-P 2 )(l-/')K- Rearranging terms in the expression of <32, gives the following line equation: <Mi-p 2 )(*i-A)+4a(* 2 -A) = o, which is exactly the equation of Line L12. Hence, when inventory is not allowed for Part Type 2, Boundary <B2 becomes Line LI 2. Corridor Windows: To determine the corridor window, we divide expressions (4.32a), (4.32b) and (4.32c) by e$ then take the limit when c* -> co and get linWc 2 + =/5T = 0; 110 ljni^/c 2 + - fi « ((1 - fl,)/2(l - /fK<t which implies that ffiw 22 + (%u> 22 + ^ = /£ > 0. Hence, the cost of reaching the cyclic schedule from Point B is always less than that of reaching the cyclic schedule from Point D. Therefore, the corridor window is the whole corridor wall given by the equation of Line L12. Figure 4.19 illustrates this policy. <B1 x 2 C x x \\ D \ 2 1 bY \ " A ^ V \ Figure 4.19: Case where c 2 — > oo. 4.6.5.4 No inventory for either part typ e In this case, we have c, + -> oo and c 2 -> oo. In a similar fashion as in the previous case, we can show that, <B1: jc, = and (B2: x 2 = 0. This means that all Ill activities take place in the third quadrant of the x-space. Also, it is not difficult to show that when c* -> °o and c* -> <», Wl -> C" = ( SA- (1-A) (1-P) Sd,,0 and IV2 -» A" = { 0,S 2 d, - (1-P) <&*, implying that the boundaries <Bi and <22 have no windows. Therefore, the optimal surplus trajectory will keep bouncing back and forth between the two boundaries and never reach the cyclic schedule in finite time. In other words, the cyclic schedule can only be reached asymptotically when Cj + -> oo and c* -» oo (Figure 4.20 illustrates this case). <B1 C x 2 *1 \ \ A <B2 Figure 4.20: Case where c£ -> oo, c* -> oo. 112 4.7 Summary In this chapter, we studied the transient behavior of a deterministic one- machine two-product manufacturing system with setup times and no setup costs. We divided the surplus/backlog space into two major Regions 91° and 9* u . For initial surplus levels in Region 9T, the optimal solution was obtained by inspection. For initial surplus levels in Region 9t u , two approaches were adopted: an algorithmic approach and an analytical one. For the algorithmic approach, a novel algorithm, to obtain the optimal state trajectory, was developed. The algorithm involves only a line search procedure, which makes it very fast. The algorithmic approach is useful for extending this work to the multi-part-type case. The analytical approach consists of the explicit expression of the optimal switching policy in Region 9T. The analytical approach is useful for studying special cases and determining the pertinent parameters of the system. Besides, it is the perfect real-time scheduling tool. The complete optimal transient solution is a feedback control policy. CHAPTER 5 OPTIMAL TRANSIENT SOLUTION OF A TWO-PART-TYPE MANUFACTURING SYSTEM WITH SETUP TIMES AND SETUP COSTS DETERMINISTIC SYSTEM 5.1 Introduction In this chapter, we extend the results of the previous chapter to the case where setup costs are nonnegligible anymore. In Section 5.2, we determine the optimal location in x-space of what we call the extracted cyclic schedule. In Section 5.3, we determine the optimal transient solution in Region 91°. In Section 5.4, we determine the optimal transient solution in Region 91". In Section 5.5, we study the case where setup times are negligible. We conclude this chapter with a summary of results. 5.2 Optimal Location of the Extracted Cyclic Schedule In Chapter 4, we studied the transient behavior of a system with zero setup costs. For a heavily loaded system (i.e., p is close to one), the condition p + p x -l>0 holds and hence, r\ and r* are equal to zero (see Theorem 3.2 of Chapter 3). In this case, there is no production at the demand rates for both part types over the cyclic schedule (Figure 3.3a of Chapter 3). In order to apply the results to the more general case of nonnegligible setup times and costs, we need to extract a schedule from the cyclic schedule of Figure 3.3d of Chapter 3 (which represents the general case). The extracted schedule has exactly the solution structure shown in Figure 3.3a. Figure 5.1 shows how to extract this schedule. 113 114 Figure 5.1: The extracted schedule. Now, we determine the location of the extracted schedule. The latter is characterized by the coordinates of the points A, B, C, and D in x-space (Figure 5.1). Let A = A) ,B = , C = l l LandD = As in Chapter 4, we denote by LI, the line containing Points D and P3 and by L2, the line containing Points B and P7 (Figure 5.1). LI and L2 are given as follows: Line Ll:d 2 p 1 (x 1 -P 31 ) + d 1 [l-p 1 )(x 2 -P 32 ) = 0; (5.1) Line L2: d^(l- p 2 )(x l -P 71 ) + d 1 p 2 (x 2 -P 72 ) = 0. (5.2) Then, we have 115 B = A- , D = C , ,4 eLl, B eL2, C eL2, and D eLl. Hence, the coordinates of Point A and Point C can be determined by solving the following system of linear equations, respectively: (d 2Pl (A-P 3l )+d J {i-p 1 )(A-P 32 ) = o, [d^(l-p 2 )(A-S 2 d 1 -P 71 ) + d^ 2 (A 2 -SA-P 72 ) = 0; ld 2 p 1 (C l -S 1 d l -P 31 ) + d 1 (l-p 1 )(C 2 -SA-P 32 ) = > [d i (l-p 2 )(C 1 -P 71 ) + d 1 p 2 (C 2 -P 72 ) = 0. The coordinates of the extracted schedule are then given by S 1 +(d 1 /d^)p 2 a 1 Q 2 - a$.-(4& - SdJ A = B = C = D = s 2 + S 2 ck+(d 2 /d l )p l a 2 (Q 1 -Sd 1 )-Q 2 p i p 2 /(l-p)) S, - SA +(4/4i)/WQi " a ^ " AKft " *0]. , s 3 +(d^/d 1 )p 1 a 2 {Q 1 -Sd l )-Q 2 p l pJ{l-p) ) <*, + ^d, + {d l /d 2 )p a o 1 (Q a - «,) - Oi A^»/(l - P)l S 2 +(d 2 /d ] ) A a 2 <? 1 - o,(l - a)(<?2 " <**.) ' s 1+ (d,/d 2 ) A o 1 ((? 2 - < 5d 2 )-Q 1 p^ 2 /(l-> 3 ) ' t S 2 - J.d, + (d,M )a«A - a 2 (l - A)(Qa -**»), (5.3) (5.4) (5.5) (5.6) where s,, S, and £> are the optimal steady state parameters (Chapter 3), and (I- p t ) cc, = (1-P) (i = l,2). Notice that, if there is no production at the demand rate for both parts (i.e., *i = and r 2 = 0), then Q x = q 1 and Q 2 = q 2 . In this case, we recover the coordinates of the cyclic schedule studied in Chapter 4 which are given by A = ,s 2 + 3 A ,B = %-*£) , Cm '*i+*A andD = S 2 J S 2 -5 l d 2/ 116 5.3 Partition of the x-Space into Two Mutually Exclusive Regions To obtain the optimal solution of the transient problem, we partition the x-space into two mutually exclusive major regions 9T and 9?° as shown in Figure 5.2. 9T and 91° are defined as in Chapter 4, but this time we use the coordinates of the extracted schedule. As in Chapter 4, we carry out the transient analysis in two parts. In the first part, we determine the optimal transient solution for all initial surplus levels in region 9t°. In the second part, the optimal transient solution for initial surplus levels in Region 91" is derived. Without loss of generality, we index the parts such that Part Type 1 is the part type with the larger setup time and setup cost (i.e., 5 Y > d 2 and fc, > fcj. Figure 5.2: Regions 91" and 91°. 117 5.4 Optimal Transient Solution in Region 91° To determine the optimal solution for initial surplus levels in Region 91°, we apply the same technique used in Chapter 4. First, we partition Region 91° into three mutually exclusive regions, G, Gl and G2. Then, we further partition Region G into eight mutually exclusive regions G11,G12, G21, G22, Hll, HI 2, H21, and H22. Figures 5.3 and 5.4 show the partition of Region 91° and G respectively. These regions are defined as follows: Gil = (K,x 2 )K -S& > 0; -d^(x 1 -P 71 ) + d 1 (x 2 -P 72 ) > 0}; G12 = {(x^KK -P 71 )-d,(x 2 -P 72 ) > 0; -d 2 (x 1 -E 1 ) + d,(x 2 -E 2 )>0; d^(l-p 2 )(x 1 -A l ) + d^ 2 (x 2 -A 2 )^0}; G21 = {(x l ,x 2 )\-d 2 (x 1 -P 3i ) + d,(x 2 -P 32 )>0; d 2 (x 1 -E 1 )-d l (x 2 -E 2 )>0; d^p i (x,-C 1 ) + d l (l-p l )(x 2 -C 2 )^0}; G22 = {(x^xjx, - S& > 0; d 2 {x,-P 3l )-d,{x 2 -P 32 ) > 0}; Hll = {(x 1 ,x 2 K- < J 1 d 1 >0; d 2 (x 1 -P 71 )-d l (x 2 -P 72 )>0; -d J (x 1 - S f 11 ) + d 1 (x 2 -sf ia )^0; -d 2 {l-p 2 ){x,-A l )-d^ 2 (x 2 -A 2 )>0}; m2 = {(x 1 ,x 2 )\d i p 1 (x i -C l )-d,(l-p 1 )(x 2 -C 2 )>0; -^K-fifn) + d 1 (x 2 -g 12 )>0; -d 2 (l-p 2 )(x 1 -A 1 )-d i p 2 (x 2 -A 2 )>0r, H21 = {(x 1 ,x 2 )\d 2 (x 1 -I 1 )-d l {x 2 -I 2 )>0; d,(l- p 2 )[x x - \)-d^ 2 {x 2 - >y>0; -d i p l (x l -C 1 )-d 1 {l-p 1 )(x 2 -C 2 )>0; 118 M , Gl D *2 / / f G \ /\? x i B V \ \ G2 \L12 Figure 5.3: Partition of Region 9t°. H22 = {(x 1 ,x 2 )|x 2 -^d ! >0; -d 2 (x 1 -P 31 )-d 1 (x 2 -P 32 )>0; -d 2 p 1 (x 1 -C 1 )-d 1 (l- A )(x 2 -C 2 )>0; -^(^ - sU +d i(*2 - g») > o}; Gl = {(x 1 ,x a )\-x l +S i d i > 0; d^x, -P 71 ) + d,(l - A )(x 2 -P 72 ) > 0}; G2 = {(Jq.xj-^ +<? 2 d a > 0; d^(l - p^x, - P 31 ) + d,p 2 (x 2 - P 32 ) > 0} , where & = [g xl ,g l2 ) T is the point in x-space given by ftl = $d, and 5l2 = C 2 +(Cj - SA'AfiiW "A). and g 2 = (g 21 ,g 22 f is the point in x-space given by 021 = A +(A " W*AKQ--P») aild ^2 = 3A- £ = (£j,£ 2 f is the point in x-space given by 119 where *2 /LG1 Gil / LE / G12 / /LG2 \ Hll /G21 / L21 ^A H12 / / / / / ^<A H21 / V ^\ \ / H22 / G22 D A / \> *i \/ \ /A B V \L12 Figure 5.4: Partition of Region G. E i = (-e 2 +^ 2 + e 1 e 3 )/o 1 ; ^-(crMX«?-«:)+k/*Xi-tf)j # 3 = 2(fc, -KM4KI& - £)+te/<Ofi; ^ = -(d,/d 2 )(l- A )/ A ; 4 = AHW)(i-a)/a; 4 = -(d 1 /d 2 )(p 2 / A ); ^ = p 2 A + WaXi - aK +(1 " A)C t + A^Mfe 4 = (!-p2)/a; 120 ^=p 2 C 2 -A 2 (l-p 2 )(l-p 1 )/p 1 +(d 2 /d l )(l-p 2 )(C 1 -A)- Remark: Line LE is the boundary in x-space separating Region G12 and G21. Since the direction of Line LE is known (given by the vector (-dp-dj in x- space), it is sufficient to find the coordinates of Point E in x-space to completely define Boundary LE. To determine Point E, we note that for initial surplus levels belonging to G12RG21 (i.e., on Line LE), the cost of reaching the cyclic schedule at Point E with a setup for Part Type 2 and the cost of reaching the cyclic schedule at Point F with a setup for Part Type 1 must be equal. Mathematically, this can be written as follows: where J represents the cost of reaching Point E with a setup for either part type. Noticing that E eLlflLE and F eL2flLE, we express F v F 2 , andE 1 as a function of E 2 and then substitute in the above cost equality. The result is a quadratic in £ 2 . £, is then obtained by substituting E 2 by its expression from the quadratic solution as shown above. If £, > A,, Region G21 becomes empty (G21 = 0). This case corresponds to a very high setup cost of Part Type 1 compared to that of Part Type 2. Therefore it is more economical to setup the machine for Part Type 2 before the cyclic schedule is reached. In the case of equal setup costs for both part types, Point E coincides with Point J. In this case, we obtain the exact partition obtained in Chapter 4 for the case of no setup costs. To summarize the control actions in Region 91°, let x = (a,b) be the vector of initial surplus levels. 121 • If x eGl: 1 : Setup the machine for Part Type 1; 2: Produce Part Type 1 at the rate £/,; 3: When the surplus level of Part Type 1 becomes 0, change the production rate to d, ; 4: When the surplus level of Part Type 2 becomes x 2 =A 2 + Ai&pjdfi - pj, switch to the control actions of the cyclic schedule. • IfxeGllUHll: 1: Do not produce either part type; 2: When the surplus level of Part Type 1 becomes 8^, start a setup for Part Type 1; 3: When the setup is over, produce Part Type 1 at the demand rate d^; 4: When the surplus level of Part Type 2 becomes x 2 = Aj + AiCLtpJdJL - a), switch to the control actions of the cyclic schedule. .If jceG12UH21: 1 : Do not produce either part type; 2: When the surplus level of Part Type 2 reaches Level l 2 , immediately start a setup for Part Type 2. l 2 is given by l 2 = (1 - p 2 ){b - (d./c^a) + (daMHl - p 2 )A +P 2 A J 3: When the setup is over, switch to the cyclic schedule control actions. • IfxeG2lUH12: 1 : Do not produce either part type; 122 2: When the surplus level of Part Type 1 reaches Level Zj , immediately start a setup for Part Type 1 . Z, is given by { = A (b-(d 2 /d 1 )a) + (d 2 /d 1 ) A C 1 +(l- A )C 2 ; 3: When the setup is over, switch to the cyclic schedule control actions. • If x eG22UH22: 1 : Do not produce either part type; 2: When the surplus level of Part Type 2 becomes ^d,, start a setup for Part Type 2; 3: When the setup is over, produce Part Type 2 at the demand rate d^; 4: When the surplus level of Part Type 1 becomes x, = Cj + C 2 d i p 2 /d 2 (l-p 2 ), switch to the control actions of the cyclic schedule. • If x eG2: 1 : Setup the machine for Part Type 2; 2: When the setup is over, produce Part Type 2 at the rate U 2 ; 3: When the surplus level of Part Type 2 becomes 0, change the production rate to d^ ; 4: When the surplus level of Part Type 1 becomes Xj = Cj + Cjd^/cyi - p 2 ), switch to the control actions of the cyclic schedule. 5.5 Optimal Transient Solution in Region 9t u As has been shown in Chapter 4, the optimal transient solution in Region 5R U is characterized by the optimal setup switching policy in that region, 123 since we know that the optimal production rates in Region 9T are the maximum production rates for both part types. In addition, in Chapter 4, it has been established that the optimal switching policy is a special corridor policy with two windows. In the case where production at the demand rates is possible in the cyclic schedule, this result is still valid and in fact is optimal when applied to the extracted schedule characterized by its location in x-space given by A, B, C, and D. Figure 5.5 illustrates this corridor policy which is completely characterized by its corridor walls defined by boundaries <31 and <B2 and its corridor windows defined by the pairs of points (W1,V1) and (W2,V2). Using the technique developed in Chapter 4 (Section 4.6), we can determine the corridor boundaries and windows as follows: Figure 5.5: Illustration of the setup switching policy. 124 5.5.1 Equations of Boundaries <B1 and <B2 (Bl: Mft + M 2 x 2 +M 3 = 0, (5.7) where M x = c-j(l-p 1 )-p i p 2 c;/(l-p 1 )(l-p)-(d l /^)(pjp 2 )a 2 c; ; (5.8a) M 2 = -K/<yAA7(l - P) ~ <IP, ~ «i(l " P2K/P2 5 (5-8d) M 3 «-(M 1 C 1 +M i C a ), (5-8c) <B2: iV,^ + N 2 x 2 + N 3 =0, (5.9) where N, = -(d 2 /d,)p 1 c + J(l-p)-c:/p 1 - a,(l -p 2 )c-J Pl ; (5.10a) n 2 = c 2 -/(i- a )- W2 c 2 + /(i-a)(i-/')-KMWa)« 1 c 1 "; (5-10b) JV 3 = -(iVA+^A)- ( 5 - 10c ) 5.5.2 Corridor windows As in Chapter 4, Let Wl = (w n ,w l2 ), Vl = (v llf v 12 ), W2 = (w 21 ,w 22 ), and V2 = (v 21 ,v 22 ) in x-space. Then, Wl satisfies the following system: \a 1 wl l +a 2 w ll +a 3 = 0; [w 12 = -(MJM 2 ]w 11 -(MJM 2 ); \v n = w u +p 1 p 2 (uj 11 - C i )/(l-p) + ((^/(^)p 2 (l-p 1 )(w 12 - C 2 )/(l-p); \v 12 = C 2 +p 1 p 2 (C 2 -w l2 )/(l-p) + (d 2 /d,)p 1 (l-p 2 )(C 1 -w n )/(l-p) (5.11) (5.12) where _ [& - l)cT ~ qfoT , [f l2 - (MJM 2 f)c* - rf l2 c 2 , 2d,(l- A ) 2^(1-^) 4(1 -A) 125 [VisMis ~ M,M 3 /M 2 )c 2 - tj 12 ( v 12 - S 2 <L,)c 2 + — ; ; (5.13b) 2d 1 (l-p 1 ) (/4 - (M 3 /M 2 f - 1>)4 - ( v 12 - S 2 d,fa 2 -K (5.13c) 2^(1-^) *hi = pMKWJMJa, -pj(l-p); (5.14a) r] 12 = «JKM)A - (MJM 2 )(1 - A )}; (5.14b) 14, = a 1 (l- A )A + «^ 2 K/ d 2 )(A+(^ 3 /^ 2 )) + ^^ 2 /(l-/'); (5-15a) v 12 = -p 1 a 2 (d^/d 1 )A l -p^ 2 A 2 /(l-p)-(l-p 2 )a 1 (MjM 2 )-a 2 {d^/d l )S 1 d l ; (5.15b) y u = ^{(I-aJ-^KMH^/MJ}; (5.16a) y 12 = p i p 2 (M 1 /M 2 )/(l-p)-p l a 2 (^/d l ); (5.16b) ftj = - Pl p 2 C 1 /(l-p)-p 2 a 1 (d 1 /cL i iC 2 +(MJM 2 )); (5.17a) ft 2 = A«aKM)Ci + (1 " A)«A + AA (M 3 /M 2 )/(l " P); (5- 17b) W2 satisfies the following system: \b<wl 2 +b 2 w 22 +b 3 =0; [w.^-iNjNJw^-iNjNJ; [v 2l = A +P&(A l -w 2l )l{l-p) + (d l l<L i )p 2 (l-p l )(A^-w 22 )l(l-p\, K 2 = *" 22 + AA^aa ~ A)/( 1 -P) + (^K)PA 1 -P2)( iu 2i ' A)/(l ~P) (5.18) (5.19) where b _ (^i ~ (NjNjft - rf 21 c[ | \y\ 2 -\)c 2 - if 22 c + 2 t (520a) 2d,(l- A ) 24,(1-^) ' *> 2 = (r 2 i^ 2 i - N 2 NjN?)c? - T} 21 (v 21 - ^Iq d,(l-A) l (r 22 (A> 2 - <W + 3AK - ^ 22 v 2 2 c 2 (5.20b) 126 *> 3 = 2d^(l- A ) + -Ki 5.20c) 24,(1 -A) 77 21 = a 1 {(d l /d I ^ 2 -(JV a /JV 1 )(l-A)}; ( 5 - 21a ) ^ 22 = Ao a (^M)(^M)-AA/(i-p); ( 5 - 21b ) v 21 = -/v^qAl-^-Aaaf^MA-a-A)^^/^)-^^/^)^; (5.22a) v 22 = Pfrfc/dJIfc + (NJNJ) + (1 - A )«A + tf 2 d^/(l - p); (5.22b) y»i = ^ 2 (^ 2 /^ 1 )/(i-^)-p 2 « 1 ( d iM); ( 5 - 23a > X 22 = a,{(l-A)-AKMKW)}; ( 5 - 23b ) / u 21 = (1 - aKA + A«ite/4M +aa WJW-A; ( 5 - 24a ) ^ = -p l a 2 [d l /d 1 iA 1 +iNjN 1 ))-p l p a A a /[l-p). (5.24b) Notice that in (5.11), w„ must be a real root of the quadratic such that u> n < C 1 . Similarly, w 22 in (5.18) must be a real root such that w 22 < A,- If either quadratic in (5.11) or (5.18) has no real root, the corridor window does not exist. Which means that, it is always optimal to reach the extracted schedule from the other side (i.e., the side with the corridor window). The optimal trajectory emanating from a point X=(a,b) in Region 5R U is obtained by using Procedure 1 (Section 4.6). Compared to the case of zero setup costs, the only change is in the expressions of a 3 and b 3 given above by (5.13c) and (5.20c) respectively, where -fcj was added to a 3 and -fc, was added to b 3 originally obtained in Chapter 4. This result is intuitive, since in the case of nonzero setup costs, we have a wider corridor windows, which means that if the setup costs are high then it might not be optimal to add another setup change and therefore reach the extracted schedule producing the same part type that the machine is set up 127 for, avoiding a rather expensive setup. Also, note that the setup costs do not appear explicitly in the equations of Boundaries <B1 and <2$2, but the latter depend on the setup costs implicitly through the coordinates of Points A and C of the extracted schedule. 5.5.3 Numerical example Consider the following system: Part Tupe 1: c, + = $2/unit/day, cf = $10/unit/day, U 1 = 8 unit/day, d\ = 2.5 unit/day, S 1 = 1 day, and fc, = $100; Part Tupe 2: c* = $2/unit/day, c~ = $10/unit/day, U 2 = 10 unit/day, ti\ = 3 unit/day, S 2 = 1 day, and fc, = $120. Figure 5.6 shows the optimal location and shape of the cyclic schedule as well as the optimal trajectory emanating at Point (-20,-20) ( i.e., we start 20 units short of both part types) in Region SR U , as well as Boundaries <B1 and <B2 and the two pairs (VH,V1) and (W2,V2) defining the two corridor windows. The optimal production and setup planning for this system can be described as follows (after rounding to integers) : Set up the machine and produce Part Type 2 at the rate of 10 units/day until its surplus level reaches 14 units (i.e., when the surplus trajectory hits Boundary <M\\ set up the machine and produce Part Type 1 at the rate of 8 units/day until its surplus level reaches 10 units (i.e., when the surplus trajectory hits Boundary <82), set up the machine and produce Part Type 2 at the rate of 10 units/day until its surplus level reaches 15 units (i.e., when the surplus trajectory enters Window (W1,V1)); set up the machine and produce Part Type 1 at the rate of 8 units/ day until its surplus 128 level reaches 16 units (which corresponds to Point P3 on the cyclic schedule); switch to the control actions of the cyclic schedule given as follows: Set up the machine and produce Part Type 2 at the rate of 10 units/day until its surplus level reaches 0, continue producing Part Type 2 at the demand rate of 3 units /day until the surplus level of Part Type 1 drops to 5 units; increase the production rate of Part Type 2 to 10 units/day until its surplus level reaches 16 units; set up the machine and produce Part Type 1 at the rate of 8 units/day until its surplus level reaches 0; continue producing Part Type 1 at the demand rate of 2.5 units /day until the surplus level of Part Type 2 drops to 8 units; increase the production rate of Part Type 1 to 8 units/ day until its surplus level reaches 16 units, which is the point we started the cyclic schedule at. 5.6 Case of Zero Setup Times and Nonzero Setup Costs In this section, we consider the case where setup times are very small or are order of magnitude smaller than the other parameters of the system (especially setup costs). This kind of situation may arise when the setups can be accomplished relatively quickly, but the equipment, the material (cleaning products for instance), and the labor used for setups are very expensive. 5.6.1 Optimal location of the extracted cyclic schedule Theorem 3.3 of Chapter 3 states that it is not possible to have simultaneously no production at the demand rate for both part types, in the case of zero setup times. Therefore, the general structure of the cyclic schedule is as shown in Figure 5.7. The coordinates of the cyclic schedule in this case are given as follows: 129 -30 -20 -10 10 Figure 5.6: Example illustration 20 Pl = S 2 +« 1 (d s /d L )/3 L /(l-A), P2 = S 2 +s 1 (d s /d 1 )p l /(l - a) - ^ P3 = P4 = v s 2y P5 130 Figure 5.7: General structure of the cyclic schedule in the case of no setup times. P6 = '$ +s 2 (d 1 /<^)pj(i-p 2 )- tA P7 = P8 = V S 2j where s,, S,., Q, and x i (i=l,2) are the optimal steady state parameters. The extracted schedule in this case reduces to a single point in x-space (Point J in Figure 5.7). To see this, let S x = 5 2 = in the coordinates of the extracted schedule calculated for the nonzero setup times case (Equations (5.3)-(5.6)). This gives A = B = ( S, + {d l ld 2 )p 2 a l Q 2 - (1 -p 2 )a 1 Q 1 ^ s 2 + {d 2 /d 1 )p 1 a 2 Q 1 - p i p 2 Q 2 /(l - p) 131 C = D = s i + (diK )whQ a - aaAA 1 - p) -n\\ S 2 + {d^/d, )p 1 a 2 Q 1 - (1 - A)a 2 <? 2 j Substituting s ; = S { - Q in the coordinates of A and B, gives Points C and D. Therefore, ( S, + (d 1 /d 2 )p 2 a 1 Q 2 -(l-p 2 )a 1 Q 1 ) A=B=C=D=I= ,s 2 + {d 2 /d l )p 1 a 2 Q 1 - p i p 2 Q 2 /(l - p) J 5.6.2 Optimal transient solution As for the case of nonzero setup times, we divide the x-space into two mutually exclusive major regions. Region 91" and Region 9T defined as follows: SR-.flx,,^) I d 2 (\-p 2 ){x i -I y ) + d i p 2 (x 2 -I 2 )<0; d 2 p 1 {x 1 -I 1 ) + d 1 (l-p 1 )(x 2 -I 2 )<0}; *° = {(x 1 ,x 2 )\ d 2 (l-p 2 )(x 1 -I 1 ) + d 1 p 2 (x 2 -I 2 )>0; d jA (x 1 -J 1 ) + d l (l- A )(x 2 -J 2 );>0}. This partition of the x-space is shown in Figure 5.8. In this case, Lines L12 and L21 coincide with Lines L2 and LI, respectively. Lines LI and L2 are defined as follows: Line LI: d^x, -I 1 )+d l Q.-fi l ){x t -I t ) = 0; Line L2: 0,(1- p 2 )(x l -I 1 ) + d i p 2 (x 2 -I 2 ) = 0. In the following, we study the optimal transient solution for initial surplus levels in Region 91". The optimal transient solution for initial surplus levels in Region 91° can be obtained by inspection in a similar way as in the case of nonzero setup times. 132 SLI P7=P8 P2 *2 PI 9T \ P5 Xy P6\ P3=P4 \L2 Figure 5.8: Partition of x-space in the case of zero setup times Recall that the optimal switching policy in Region 9T for the case of nonzero setup times is a special corridor policy characterized by boundaries <B1 and <B2, and the two windows defined by the pairs (W1,V1) and (W2,V2). In the case of zero setup times, to obtain the optimal policy, all we need to do is to let <?! = and 8 2 = in the equations of Boundaries <B1 and <B2 as well as in the coordinates of Points Wl, W2, VI and V2. Equations of Boundaries <B1 and <B2 (Bl: MjXj + M 2 x 2 + M 3 = 0; (B2: N 1 x 1 +N 2 x 2 + N 3 =0, where M x , M 2 , N x , and N 2 are as defined by (5.8a), (5.8b), (5.10a) and (5.10b) respectively, and M 3 = -(MJ, +M 2 I 2 ), and N 3 = -{N& +N 2 I 2 ). 133 Corridor windows Wl and VI are given by ia 1 w^+a 2 w n +a 3 = 0; \w 12 = -(MjM 2 )w n -(MjM 2 ); \v u = w n + Pl p 2 (w u - CJ/fl - p) + {a\/a\)p 2 (l - A )(u/ 12 - C 2 )/(l - p); \v 12 = C 2 + p i p 2 (C 2 -w 12 )/(l-p) + (a\/a\)p 1 (l-p 2 )(C 1 -w ll )/(l-p); W2 and V2 are given by ib 1 w 2 l 2 +b 2 uj 22 +b 3 =0; iv 21 = A + AA(A " w 21 )/(l - />) + (d,/d,)/o,(l - aMA, - w 22 )/(l-p); \v 22 = w 22 + aa (">* " As )/(l - P) + (4iM )A (1 " A )Ki " A )/(l - A* where, Oj, a 2 , a 3 , b x , b 2 , and b 3 are obtained by plugging £,=0 and £ 2 = in (5.13a), (5.13b), (5.13c), (5.20a), (5.20b) and (5.20c) respectively. Numerical Example We consider the same example studied in Section 5.5 and let the setup times be zero and double the setup cost of each part type. The reason for increasing the setup costs for both part types can be thought of as using more sophisticated and costlier techniques, tools or labor to reduce the setup times to a nonsignificant level compared to the other parameters of the system. Part Tupe 1: c; = $2/unit/day, c," = $10/unit/day, £7, = 8 unit/day, d\ = 2.5 unit/day, <?, = 0, and fc, = $200; Part Time 2: c 2 m $2/unit/day, c 2 = $10/unit/day, U 2 = 10 unit/day, d\ = 3 unit/day, S 2 = 0, and /^ = $240. 134 The optimal trajectory and the optimal shape and location of the cyclic schedule are shown in Figure 5.9. Figure 5.9: Example illustration 5.7 Summary In this chapter, we extended the results of Chapter 4 to the more general case of nonzero setup times and nonzero setup costs. We showed that the optimal trajectory reaches the extracted schedule first, then follows the a cyclic schedule thereafter. Also, we studied the case where setup times can be 135 neglected compared to the remaining parameters of the system. In this case, we showed that the extracted schedule reduces to a single point in x-space, and obtained the optimal transient solution analytically as a special case of the more general case of nonzero setup times and costs. CHAPTER 6 OPTIMAL AND NEAR OPTIMAL CONTROL OF A TWO-PART-TYPE STOCHASTIC MANUFACTURING SYSTEM WITH NONRESUMABLE SETUPS 6.1 Introduction In this chapter, we consider the same system studied in the previous chapters with the following additional assumptions: The manufacturing system is unreliable, the setup times are now random and nonresumable, the planning horizon is infinite, and the optimization criterion is the total discounted cost over the planning horizon. The purpose is (ij to propose a stochastic control model for the simultaneous planning of production and setups in the unreliable manufacturing system; (ii) to develop an efficient technique for the computation of the optimal control policy; and (in) to provide heuristics that are efficient and simpler to implement in the real world. 6.2 Mathematical Formulation Consider a manufacturing system with a single machine capable of producing two kinds of parts. The system should satisfy a constant demand rate for each part type. The machine is unreliable and is subject to random breakdowns and repairs. We assume that failures and repairs are time dependent (see Gershwin, 1993). The time to failure and time to repair are modeled as exponentially distributed random variables. The machine is not perfectly flexible in the sense that changeover or setup times between production of different part types are significant. That is, 136 137 the machine incurs a nonzero setup time when switching from one product to the other. The setup times are assumed random with exponential distributions. Indeed, In many manufacturing systems, a setup change requires the calibration of the machine and some parts that fail inspection may be initially produced. Thus, requiring a random amount of time. Also, the assumption of a random setup time is very convenient mathematically as will be shown later. The machine could be either working or under repair. If the machine breaks down, it is immediately brought to repair. During a repair the machine can neither produce parts nor change setups. The machine may breakdown during a setup. This assumption is more realistic than that of no failure of the machine during a setup (see Hu et al., 1995), especially that we consider random setup times. Usually, setups are assumed resumable. That is, in the case of an unreliable system, after a repair, the setup resumes from the point it was left in, just before the failure occurred. This assumption is not always true, especially if the repair involves disassembling parts of the machine or when the manufactured products are chemical fluids for example. In this dissertation, we assume that the setups are nonresumable. 6.2.1 System dynamics The system considered has a hybrid state comprising both a continuous and a discrete component. The continuous state equations Let x,(t) be the cumulative production surplus of Part Type i (i-1,2) at time t. A positive value of x,.(t) represents inventory, while a negative value represents backlog. 138 Let u.(t) be the controlled production rate of the machine producing Type i parts at time t. Let d i (i=l,2) be a given demand rate for Part Type i. The rate of change of the surplus of Part Type i is given by dx(t) dt = u(t) - d, (6.1) where x(t) = ( Xl (t),x 2 (tj), u(t) = (u,(t),u,(t)), and d = (c^d,). The discrete event stochastic process The discrete part of the state variable describes the system's operational mode. Let M = {0,1,2,3,4,5}. The operational mode of the machine is described by the random variable £(f) with values in M, which is thus the value at time t of a controlled jump process with six possible states. The meaning of each possible state is given in Table 6.1. Let A. be the transition rate from state t e M to state jeMat time t. /L = lim — Pr{<f(f + dt) = j\ #t) = i}, for {i,j)eMxM ••«•«■ dt (6.2) Let Q \tr the infinitesimal generator of the stochastic process {£(i), t £ 0} Q is then given by -r r P -P- ~ K* -Ks *14 ^15 P -P - ^25 ^25 P -p-^34 ^34 P «1 "P-Sj P S 2 -p- (6.3) 139 Table 6.1: Operational mode description ffl Meaning Machine under repair 1 Machine up but not set up 2 Machine changing setup to Part Type 1 3 Machine changing setup to Part Type 2 4 Machine set up for Part Typ e 1 5 Machine set up for Part Type 2 where p and r are the failure and the repair rates of the machine, respectively. Sj and s 2 are the rates of changing setup from Part Type 1 to Part Type 2 and from Part Type 2 to Part Type 1, respectively, p, r, s lf and s 2 are given constants. The transition rates X^, X 2S , A 14 , and A 1S are all control variables. These are modeled as impulsive controls which trigger an immediate jump of the controlled process (see Haurie and L'Ecuyer (1986)). That is, they can take values only on the set {0,oo}. We model these controls as actions through which we decide whether to initiate a new setup for a part type or to continue with the current setup. Let A be the set of actions that the controller can choose from, which is defined by A = {a ,a lt a 2 }. These actions are described in Table 6.2. A state dependent constraint on the action set is defined by a subset T of MxA, called the set of admissible state-action couples, such that AW)) = {«(*) ** I mMt)) er} V %t) eM. (6.4) Table 6.2: Description of the controller action Action Meaning °0 Do nothing ^1 Initiate a setup change for Part Type 1 a 2 Initiate a setup change for Part Type 2 140 Analytically, the set r is given by T = {(0,a ),(l,a ),(l,a 1 ),(l,a 2 ),(2,a ),(2,a 2 ),(3,a ),(3,a 1 ),(4,a ),(5,a )} (6.5) The complete state of the system is given by the vector (x(f),£(t)) e SR 2 x M. The vector v(t) = (u(t),a(t)) defines the complete control of the system. Where, u(t) = (ujt^u^t)) specifies the production rates of the two part types at time t, and a(t) specifies the setup decision at time t. 6.2.2 Control constraints Let V i denote the maximum production rate of the machine when it is producing Type i parts. Then, we have the following two constraints: 0<u l <U,I [mm2] ; (6.6a) O^u^Ly^,, (6.6b) where [1, if £(t) = i, ieM {WN} [0, otherwise. The complete control space is defined by the set S(Kf)) = {(u(t),a(t)) I < H < t/,./ U(tH+1) , i=l, 2; and a{t) e^t)}, (6.7) which is a compact set. 6.2.3 Penalty function We assume that a cost rate of % dollars per unit time (^1,2) is incurred when the machine is undergoing a setup from Part Type j to Part Type i {jtij . The instantaneous cost which penalizes the production for being ahead of (i.e., jc ; >0) or being behind (i.e., x^O) the demand is given by £" (c/xf + cpefj. 141 Where, c* and c," are the per unit instantaneous inventory holding and backlog costs respectively, and x*(t) = max{x i (t),0} and x7(t) = max{-x i (t),0}. The total instantaneous cost is then given by 9(x, &v) = X';' (« + c i x i) + Xj{t-2) + Mf-a) • < 6 - 8) 6.2.4 The stochastic control problem Admissible policies Let (3 t :t > 0} be an increasing class of sigma-fields describing the history of the process {(*(<), £(t)):i > o}, and let (£2,3) be a measurable space. Then, a sample realization a> eQ. corresponds to a trajectory of {x(t):ti:0} which is continuous and to a sequence of values of {£(*): t £ 0}. The set II of admissible control policies is a family of 3, - adapted processes with values in E(Q, such that the mapping n : tt 2 x M -> E(£) n{x,%) h-» v = (u,a) (6.9) is piecewise continuously differentiable (e.g., sufficiently smooth). Thus, with each control policy, n ell, is associated a probability measure P n on (£2,3) such that the process {(x(t),^{t)):t > 0} is well defined. An admissible control policy is a feedback control that specifies the control actions when the system is in a given state. The feedback control determines the production levels and the setup state of the machine as a function of the surplus level and the status of the machine. Cost functional The objective is to find a control policy n* ell, corresponding to a production flow rate control u* = (u^,^) and a setup control a*, which 142 minimizes for each initial condition (x ,| ) the following expected discounted cost JMo^o) = E x {[e->*g(x[t),fflMm\x[0) - ^O) = &}, (6.10) where the minimization is over all functions ^ x( r), £( r)) = (u(r),a(r)), such that x(r), £(r), u(r) and a(r) satisfy the system dynamics and (u(r),a(r))e3(^) for r e[0,+oo). Here, p is a given positive discount rate. 6.2.5 Optimality necessary and sufficient conditions This optimal control problem belongs to the class of problems considered by Rishel (1975). It is to control a stochastic process with jump Markov disturbances. A necessary and sufficient condition characterizing the optimal control policy, under appropriate differentiability assumptions of the control (Rishel (1975)), are described by the following set of partial differential equations (Known as the Hamilton- Jacobi-Bellman (HJB) optimality equations): pJ(x,Q = rnm\g( X ,Z,v) + f j ^^x 1 + Z%A v y j ( x >?)}> V ^ eM ( 6 - n > B-B K>[ ,--i dx { ?eM J where J(x,£) denotes the optimal cost value. q^(v) is the transition rate from state | to state £ , which depends on the control v . Here, x t = u, - d, denotes the total derivative of x, with respect to time. Notice that the production rates u t (i = 1,2) should be chosen so as to minimize the rate at which the cost increases with respect to x i . The setup decisions are determined independently from the u/s. Based on Rishel's (1975) formalism, the optimal control formulation has been applied to a number of flexible manufacturing systems (Tsitsiklis (1982), Akella 8s Kumar (1986), and Kimemia and Gershwin (1983)). Although there is 143 no available solution technique to solve this control problem analytically, the previous work has revealed some characteristics of the optimal solution. A general finding is that the optimal solution is a switching policy. That is, the control policy is determined by a partition of the state space, and control actions take place when the system state moves from one region to another. In the next section, we adopt a numerical approach to obtain the optimal solution of the control problem. We use an approximation technique due to Kushner and Dupuis (1992) to discretize the HJB partial differential equations. This technique leads to an infinite horizon stochastic dynamic programming algorithm which we solve using the successive approximation technique. Then, based on the numerical solution, we observe the structure of the optimal solution and construct a heuristic for real-time scheduling. 6.3 Numerical Approach 6.3.1 Markov chain approximation method In this section, we adopt Kushner and Dupuis's (1992) method to construct an approximating discrete time discrete state Markov chain to the original continuous stochastic control problem, we have deliberately tried to keep the notation as close as possible to that used by Kushner and Dupuis. The general idea of Kushner and Dupuis's approach is to approximate the value function J[x,%\ by J h [x,%\ via a finite difference method. Basically, the partial derivatives dJ\x i ,^\jdx i are approximated by the following one sided difference approximations: 144 where h t e 9t + represents the finite difference interval of the variable x i . e, e <R 2 is the unit vector in the i* coordinate direction. Notice here that the direction (forward or backward) of the finite difference approximation (6.12) is the direction of the production surplus velocity component x i . Let b i (x) = x i (t) and define the following quantities: b*(x) = maxfbJjc^O) and b7(x) = maxf-b^x^O). Then, aT{x,$ _ J fc (x+frq,g-j fc (x,fl y dx + J h M -J h (x- hi e^ b . (x) (6.13) Substituting in (6.1 1) and using the fact that q^v) = -V q^^K we obtain f J [x,£] = muv !=*(b;tx)+b:ix)) 1-2 1-1 h. h,. + £q>y h (*,<r) ; <?*•? Let Q h (x,$,v)= X^M + Zflfr.Ml/M- Noticing that &+(*)+ b^H^N' dividing <•#* i-i through by £) h (x,£,v) in the previous equation and rearranging terms, gives > i + Q h (x,4,v) J h (x,3 = min. •>*M[Q n (x,Z,v) b > (x) ./(x + V 1 ^) + -* r /(x-M^) Gr[x,&vy\ g h (x,^,y)h, Define the grid S h where each component of the surplus vector x, can only take values that are positive or negative multiples of h. . That is, 145 S h ={x\x = 2A«A : ", = 0,±1,±2,..., i = 1,2}, Let p h (x,$?,v) = * ?{V) and let Q (x,§,v) p h (x,x',£,v) = < Then, we obtain bf[x) <?(x,4,v]h l if x' = x±h i e i for all other*' eS" J h (x,£) = min- fe5<<0 g(x,4,v) Q h (x,£,v) 1 + - <?*(*,£"). 1 1+ > - x'eS* ?#* (6.14) Equations (6.14) have a nice interpretation. Indeed, they represent the dynamic programming equations of a discounted cost infinite horizon discrete time discrete state Markov Decision Process. Where, p h (x,x',£,v) and p h (x,%,%,v) are the transition probabilities (since they are nonnegative and sum to 1 over (x',<?) 6^^^ for each [x,$ eJl^S* , and I /(l+p/Q h (x,£,v)) <1 is the discount factor. In Kushner and Dupuis (1992), it is shown that this approximating Markov chain satisfies the contraction properties which guarantee the existence of a solution which can be obtained by classical iterative methods in policy or value space. Also, it is shown that with appropriate boundary conditions hmJ h (x,£) = J(x,$. h->0 146 6.3.2 Computational algorithm: The value iteration method For any function J h :Yl^ M S h -> 9t, n W itM S h ~* 3 (^ with ?rer1 ' for ^ (x,^) ef]^^, we denote: ^•^rnS V Q h (x,z,4x,m + Y,P h (x,x\ZMx,®V h {x\t) + Zp h (x,S,?Mx,®V h (x,m; (6.15) Note that T ff (j h )(.,.) is a function denned on JJ^ S h , and T„ may be viewed as a mapping that transforms a function J h on IT^S* into another function T„\J h ) on T T S h . The mapping T x plays an important theoretical role and provides a convenient shorthand notation in expressions that would be too complicated to write otherwise, from (6.15) it can be seen that T„(j h )(x,$ represents the optimal cost corresponding to policy {n,n,...} for the one-stage problem with initial state (x,§, stage cost g{x,^4x,^)/(p+Q h {x,^4x,^)), and terminal cost function :- e S k +%p h (x,Z,Z,4x>®V h {x>t) Let r(j h ){x,® = mm{Tp h )(x,^\ Algorithm Given a finite difference h e 9i 2 + : Initialization step: Choose e e5R + ; 2>"(x,x',£*(x,3)J h (x',<f) l{p+<?{x,$,n{x,Q)). 147 Let fc:=l, and {j h f(x,^:= g{x,l,(0,a ))/{l - e -») , V(x,i\ e]J^S h ; Go to main step. Main step: Let (j h ) k -\x,®:= {j h T(x,$, V(x,$ e]J^S h ; determine the policy n* such that: (j h ) k (x,®:=Ts{j h ) k -\x,3:=T\j h ) k -\x,3,V(x,®e]] 4eM S h ; Let *' V( ,^^^ )fc(X ^- (jhr(X '^ ; Go to Convergence test step. Convergence test step: lf\c k -c k \>£, Let k:= k + 1 and go to main step; else ri = 7?, stop. 6.3.3 Implementation considerations Since Grid S h is finite, we have to impose certain boundary conditions on the states of the system that fall outside Grid S h . These boundary conditions constitute an additional approximation to the original problem. However, as will be shown from the subsequent numerical examples, they do not affect the optimal solution. For each £ e M, whenever a component x i of x is on the boundary, the boundary conditions consist of approximating the cost of the state to which the system is going to jump, which is outside Grid S h , by a linear interpolation between the current state and another state which is inside Grid S h . To better understand these boundary conditions, consider the following situation: Assume that x 1 is on the boundary of S h and that x 1 +h 1 is 148 outside S h . Hence, x, -ft, is inside S h . A linear interpolation, using J h [x l ,x a ,^ and J h [x Y -\,x 2 ,%\, gives the following approximation for J h (x, + h 1 ,x 2 ,$: J h {x 1 +h 1 ,x 2 ,® = 2J h (x 1 ,x 2 ,$-J h (x 1 -f h ,x 2 ,$. In the following section, we illustrate this numerical approach by means of two examples. 6.4 Optimal Policy To characterize the optimal policy, we apply the above algorithm to two cases which data are given in Table 6.3. For both cases, the discount rate p=0.1. The grid S h is defined as follows: S h = {x\x = XAe.n, : h, = 0.5, n,. = -30,...,0,...20, i = 1,2}. Table 6.3: Data for cases 1 and 2. Case 1 Case 2 p 0.05 0.05 r 0.9 0.9 te.dj (0.25,0.25) (0.25,0.25) 1PM (0.67,0.67) (0.67,0.67) K»*a) (0.5,0.5) (0.5,0.5) (cr.cr) (1.0,10.0) (1.5,15.0) K.Ca) (1.0,10.0) (1.0,10.0) l.£l J %2 1 (1.0,1.0) (2.0,1.0) The results are shown in Figures 6.1. a, 6.1.b, 6.1.c, 6.1.d and 6.3 for Case 1 and Figures 6.2.a, 6.2.b, 6.2.c, 6.2.d and 6.4 for Case 2. 149 6.4.1 Optimal setup policy when the machine is up As shown in Figures 6.1. a and 6. 2. a, the optimal setup policy partition the state space (for E, = 2 or £=3) into four mutually exclusive regions. These are described below. Region I: If the machine is already set up for the production of Part Type 1, then it is optimal to continue with this setup and produce this part type at the maximum production rate. On the other hand, If the machine is set up for the production of Part Type 2, then it is optimal to setup for Part Type 1 immediately and start producing this part type at maximum production rate. Notice the difference in the shape of Region I in Cases 1 and 2. This is mainly due to the fact that Part Type 1 needs to be produced more often than Part Type 2, since its backlog, inventory, and setup costs are higher than those of Part Type 2. Region II: The policy in this region is just the opposite of that of Region I . That is, a setup is optimal only if the machine is already set up for the production of Part Type 1 . The production strategy is to produce at maximum production rate when the machine is set up for the production of the proper part type. Region III: The policy in this region is to keep the same setup of the machine and set the production rate to zero. In other words the production is stopped for either part type. Region IV: If the system is in this region, the production continues with the current setup at maximum rate. Regions I, II, and III are transient and Region IV is recurrent (i.e., the surplus vector (x,,x 2 ) will be trapped with probability one in Region IV, once the system moves into that region). To see this, consider Figure 6.5, where different possible trajectories of the production 150 surplus vector are shown to endup in Region IV no matter where the initial production surplus originated. Trajectory 1 (thick solid line) shows a surplus starting in Region II, with the machine set up for Part Type 2. The policy of Region II dictates that we should start a setup immediately and produce Part Type 1 at the maximum production rate. This will drive the surplus trajectory to Region IV. Trajectories 3 and 4 (dashed and dotted lines respectively) show a starting surplus in Region III, in which the policy dictates that we stop production (surplus levels too high) . The surplus of both parts decreases and either boundary of Region II (trajectory 3) or boundary of Region I (trajectory 4) is encountered. At this point, a chattering on the boundary occurs, since we produce and stop the machine alternately. This chattering continuous until the surplus vector enters Region IV, where it stays thereafter. Trajectory 2 has the same behavior as Trajectory 1 . 6.4.2 Optimal setup policy right after repair As shown in Figures 6.1.b and 6.2.b it is optimal to start a setup for the production of Part Type 1 (Part Type 2) if the production surplus vector is in Region I (Region II). If the production surplus vector is in Region III, then it is better to delay the setup until the production surplus vector, driven by the demand, drifts to either Region I or II. 151 a 10 5 -5 -10 -15 IV -10 X1 Figure 6.1. a: Optimal setup policy (machine is up) a 10 5 -5 -10 -15 I I 1 'y£.-£. ■".■..■',-;;;.•;;.•. Ill; v^vX0K : ; : " mil : :: :: !: :: :: :: ; IJI IP _ ; r ; : : : ; : : : r *-r :; : .:-■-: 111 ....,-,;,;:. mm HI ■:o:--::-"- -10 10 X1 Figure 6.1.b: Optimal setup policy right after repair Figure 6.1.c: Optimal production rate for Part Type 1 Figure 6.1.d: Optimal production rate for Part Type 2 Figure 6.1: Optimal policy for Case 1 152 a 10 5 -10 -15 - «-:-:-:■:■:■:■:.:•:■:■:■:■;-:•;■.■-;-;■;-;■:■>:•::->:>■■::■: .... :■:-:-:-»;■:■;•;■:•:-.■■»■.■;•;■;■:■;-:■:• -: :■ -:;■ ■: :■ •: ■:<:•:■:«:•:«•:•:■«■:■:■;■;•;•:■:•;•:.:■:■:•:•;•:■:■.■: ;;- ■: :■ .. .- w:-:v:v;v;--;v;v,;v;v:v:..-;:--::. -. .- -. .. -. . ■:■:■:•:■:•:'>:•:■-.■>;■;•;-:■:■:■:■:•;-:■:■:■:■:-:■::.-::.•::.-;:• :«•>:•;•;•:■;•;■;>:•;•:■:•;■:■:■:•:•:-:-:■:.:.:■:..■::•..: :■ ■: .-.■: . ■ , ■■;;■:.■.':■;'.■.■ ;-:■:•;■;•>.■.•: ■.■.■.-,■...■,■.■ ■.■■;.■.■■■;-.:■:-■-.'■ 111 -10 X1 10 Figure 6. 2. a: Optimal setup policy (machine is up) Figure 6.2.b: Optimal setup policy right after repair 0.5- Figure 6.2.c: Optimal production rate for Part Type 1 Figure 6.2.d: Optimal production rate for Part Type 2 Figure 6.2: Optimal policy for Case 2 153 6.4.3 Optimal production policy Figures 6.1.c, 6.1.d, 6.2. c and 6.2.d show the optimal production rates of the machine (when it is possible to produce) for Cases 1 and 2. As can be seen, the optimal policy is of the bang-bang control type. That is, the production rate is either at its lower (zero in this case) or at its maximum level, and there is no optimal production level in between. 6.4.4 Optimal cost function In Figures 6.3 and 6.4, we give, for Case 1 and Case 2 respectively, the optimal cost as a function of the surplus vector, as well as its corresponding contour lines for each setup (i.e., when the machine is set up to produce Part Type 1 and when it is set up to produce Part Type 2). Notice that the optimal cost function is convex (the level sets are convex). Also, it is very important to notice that the cost function is everywhere continuously differentiable which is the basic assumption for the HJB partial differential equations. 6.4.5 Policy structure In both cases the optimal policy is of the corridor type as proposed by Sharimia et al. (1991) (Region IV is the corridor). The sample trajectory shown in Figure 6.6 shows that the production surplus vector will converge to an average cyclic schedule. This cyclic schedule was referred to as the Limit Cycle in Sharimia et al. (1991). In the following, we formally define the two types of optimal policy structures obtained. Definition 6.1: We call an Orthogonal Corridor Policy (OCP), a policy having the setup policy shown in Figure 6.1. a (Case 1). 154 2000 1000 2000 1000 CNJ X 8 10 5 -5 -10 -15 10 5 -5 -10 -15 \ / I ' 1 t v^w. \\\vvv\ nSNSS: z s NNN $^ ■ -10 -10 X1 10 I I 1 IL v\X^\\\^_ __ XS^NN^- vvv^\ \^NS^^— ----^ NX\V \\/0\/O\^^ S^v^ \\\ 10 X1 Figure 6.3 : j(x 1 ,x 2 ,2) (upper 2 graphs) and j(x 1 ,x 2 ,3) (lower 2 graphs) for Case 1 . 155 2000- 1000 a 2000 1000 X Figure 6.4 : j{x 1 ,x 2 ,2) (upper 2 graphs) and j(x 1} x 2 ,3) (lower 2 graphs) for Case 2. 156 Definition 6.2: We call a Parallel Corridor Policy (PCP), a policy having the setup policy shown in Figure 6. 2. a (Case 2). In the next section, we give a necessary and sufficient condition for the demand to be feasible. Then, we propose a heuristic based on the results of the deterministic system studied in Chapter 5 to approximate the optimal policy of sufficiently reliable stochastic systems. Region III Region IV Region I Figure 6.5: Production surplus trajectories emanating in different regions. 157 Surplus trajectory x 2 Average cyclic schedule ^\ i, Vv\ /. B Figure 6.6: Typical surplus trajectory within the corridor. 6.5 Demand Feasibility In this section, we provide a necessary and sufficient condition for the demand to be feasible. Applying the conditions stated below allows us to distinguish between systems with enough capacity to meet the demand and those for which the demand will never be met. In the latter case, no matter which policy we use, the cumulative production surplus of the system will eventually drift to — 00. The following proposition gives a necessary and sufficient condition for the system to meet the demand. Proposition 6.1 (Feasibilitu Testj: The demand set, d = {d it cL i }, is feasible if and only if 158 U, _r Y a, p + r p+s ( > d 1 +d 2 , (*=1,2) Proof: The proof is based on a setup policy referred to as the Maximum Capacity Policy (MCP). Let T be the planning horizon and without loss of generality, we assume that the machine is operational and has not been set up at time t = . The MCP operates as follows: Start with a setup for Part Type 1 at time t = 0, then keep producing this type at maximum production rate U v until time t = T 1 . During the time interval (0, 7}), the machine is set up for Part Type 1 after each repair (since the setup is nonresumable) . At time t = T iy we switch to Part Type 2 and continue producing this type at maximum production rate U 2 , until time T = 7] +T 2 . During the time interval [T v T), the machine is set up for Part Type 2 after each repair. To maximize the production capacity of the T, d. system, under such a control policy, we choose 7} and T 2 so as — = — . T 2 d, Therefore, we have T. = — — — T and T, = — — — T . dt+d, d 1 +d 2 Note that, operating the system in this fashion, the number of setups is minimum and the machine is always producing at its maximum production rate when it is up. Therefore, it allows the maximum attainable system capacity. We want to show that under MCP, the condition above is necessary and sufficient for the demand to be feasible. Now, let N(t), t > be the number of machine failures at time t, i.e., {N(t), r > 0} is a renewal process with inter-renewal times t,. = t" +tf [i = 1,2,...). Where if, if (i ■ 1,2,...) are exponentially distributed with rates p and r respectively. Also, let P t (T) (£=1,2) be the cumulative production at time Tof Part Type i. Then we have 159 (Nft) P { (T)= Xtf-JTOA.fr,, (*=1,2). \i-0 Let A,, be the random variable denoting the total time spent to complete a setup. A. includes the random setup time and the time spent in repair when a breakdown occurs during a setup. The demand is feasible if and only if the following holds: lim^teM > di> [Uh2 ). r-»«c J 1 <=> E lim- 7-*oo N(r,) i-0 -£[jt(i;)]b[aj * T7>^, u.. 2). <» Eft"] = — and E[A,] = ~^— P P + s t ]im E[N[T i )] rT l > s, A p(p + s ( J7. , (i=l,2). <^> r— r 4 (i=l,2). Note that, when T -> », 7}-K». Using the elementary renewal Theorem of Renewal Theory, we have lim L J = — . Where // = — + -. The result follows Tr»» T { n p r immediately. ■ Remark: Intuitively, if the condition holds, then we know that there is at least one setup policy (the MCP precisely) for which the demand is feasible. On the other hand, if Proposition 6.1 is violated then we know for sure that there is 160 no other policy for which the demand will be feasible, since the MCP which allows the maximum attainable capacity, failed to satisfy the feasibility test. 6.6 Near Optimal Real-Time Policy for Sufficiently Reliable Systems In this section, we propose a heuristic, implementable in real-time, for stochastic systems that are sufficiently reliable. The heuristic is based on the results of the deterministic system studied in the previous chapters, it is referred to as the hedging corridor policy which is characterized by the two corridor edges defined by z, and z 2 as shown in Figure 6.7. Notice that, before we switch to the production of the other part type, a certain inventory is built to hedge against future shortages brought about machine failures and setups. 6.6.1 Construction of the hedging corridor policy The hedging corridor policy is constructed as follows: We notice, from the numerical optimal solution, that all trajectories converge to a certain average cyclic schedule and that there is no production at the demand rate. Hence, there must be no windows in the corridor (since the trajectories will never reach the cyclic schedule in finite time), and the constructed corridor should be positioned on the extracted schedule rather than on the cyclic schedule itself (see Chapter 5) . Hence, the hedging corridor policy is completely characterized by z 1 =A 1 andz 2 = C 2 , where \ is the abscissa of Point A, and C 2 is the ordinate of Point C defined in Chapter 5. To be able to use the formulas of A 1 and C 2 , we substitute the setup times d x and S 2 by the average setup times of the stochastic system, given by S x = — and S 2 = — . 161 Figure 6.7: Typical surplus trajectory within the hedging corridor. 6.6.2 Performance of the hedging corridor policy Study cases To test the performance of the hedging corridor policy, we simulated the system using the optimal numerical policies and then compared with the results using the heuristic. We used the cases shown in Table 6.4. Case 1 is referred to as the symmetric case, since both part types have the same characteristics. Compared to the symmetric case: Cases 2 and 3 show an in increase of the setup time of Part Type 1 . Cases 4,5, and 6 show an increase of the inventory cost rate of Part Type 1 . Cases 7, 8, and 9 show an increase of the backlog cost rate of Part Type 1 . Cases 10 and 1 1 show an increase of the setup cost rate of Part Type 1. 162 Cases 12, 13, and 14 show an increase of the production rate of Part Type 1. Cases 15, 16, and 17 show an increase of the demand rate of Part Type 1. For Cases 1 through 17, we used the following data: p= 0.1, p = 0.05, r = 0.9 and a grid S h = {x\x = £A € A : h > " ° 5 ' n . = -30,...,0,...20, i = 1,2}. To test the sensitivity of the hedging corridor policy to the failure rate of the machine, we increased p using the symmetric case (Case 1). p = 0.1 corresponds to Case 18. p = 0.15 corresponds to Case 19. p = 0.2 corresponds to Case 20. Performance measures • The total discounted cost of operating the system over the simulation horizon. • The average surplus level X; for Part Type i over the simulation horizon. The closer the average surplus to zero, the better the policy. • The production percentage, given as follows: f> = ^-xl00%, fori = 1,2, where W t and D. are the cumulative production and cumulative demand of Part Type i over the simulation horizon. This is the fraction of the satisfied demand of Part Type i expressed in percentage. The closer this measure to 100%, the better the policy. • The aggregate production percentage given as follows: P = t£* xl00%- • The total average surplus X, which is the sum of the individual part type average surpluses. 163 Table 6.4: Simulation cases Cases c r < c o < Xi X 2 ^ u 2 4 d, s i S 2 1 15 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 1 2 15 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 0.75 3 15 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 0.5 4 15 2 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 5 15 2.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 6 15 3 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 7 20 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 8 25 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 9 30 1.5 15 1.5 0.5 0.5 1.13 1.13 0.32 0.32 10 15 1.5 15 1.5 0.75 0.5 1.13 1.13 0.32 0.32 11 15 1.5 15 1.5 1 0.5 1.13 1.13 0.32 0.32 12 15 1.5 15 1.5 0.5 0.5 1.43 1.13 0.32 0.32 13 15 1.5 15 1.5 0.5 0.5 1.83 1.13 0.32 0.32 14 15 1 5 15 1.5 0.5 0.5 2.26 1.13 0.32 0.32 15 15 i 1 5. 15 1.5 0.5 0.5 1.13 1.13 0.4 0.32 16 15 15 15 1.5 0.5 0.5 1.13 1.13 0.48 0.32 17 15 1.5 15 1.5 0.5 0.5 1.13 1.13 0.64 0.32 1 Simulation A continuous simulation model was implemented. For each of the 20 cases above, the simulation was replicated 10 times using different random variate streams (to reduce biases). The simulation horizon (for one replication) was fixed to 10,000 time units. The initial surplus values were equal to zero in all the cases. 164 Simulation results Tables 6.5 and 6.6 show the simulation results using the optimal policies and the hedging corridor policy, respectively. Besides the performance measures defined above, Tables 6.5 and 6.6 show the structure of the optimal policy and the result of the feasibility test for each case. Table 6.6 shows that the hedging corridor policy performs very well compared to the optimal policy. Notice that for all the cases the demand is met at 100% and that the average production surplus levels are practically 0. The maximum percentage gap between the total discounted cost of the optimal policies and the heuristic ones is 23.9% and occurs with Case 17 for which the demand of Part Type 1 is double the demand of Part type 2. We might conclude that for sufficiently reliable systems producing part types that are similar in characteristics (which is usually the case in practice), the hedging corridor policy performs very well. The latter also performs well with sufficiently reliable and moderately loaded systems. When the system becomes heavily loaded or frequently unavailable the heuristic performs poorly. Notice that for some cases, the total discounted cost of the heuristic policy outperformed that of the optimal policy. This is mainly due to the continuous implementation of the simulation. Indeed, when simulating the system using the hedging corridor policy, we calculate exactly the quantity of a part type to be produced before we switch to the other. While, when simulating the system using the optimal policy, we calculate the quantity to be produced during one time unit and then check whether we hit the switching boundary or not. This might cause an overshooting of the boundary resulting in higher cost (especially, when the production rates are high compared to the demand rates) . 165 r 'able 6.5: Simulation results i using the numerical optimal policy Cases Demand Feasibility Policy Structure Total cost Xj X 2 X p* p 2 p 1 Feasible OCP 100.7 0.0 0.0 0.0 100.0 100.0 100.0 2 Feasible OCP 114.3 0.0 0.0 0.0 100.0 100.0 100.0 3 Feasible OCP 149.9 0.0 0.0 0.0 100.0 100.0 100.0 4 Feasible OCP 103.2 0.0 0.0 0.0 100.0 100.0 100.0 5 Feasible OCP 103.6 0.0 0.0 0.0 100.0 100.0 100.0 6 Feasible OCP 105.1 0.0 0.0 0.0 100.0 100.0 100.0 7 Feasible PCP 107.4 0.0 0.0 0.0 100.0 100.0 100.0 8 Feasible PCP 111.8 0.0 0.0 0.0 100.0 100.0 100.0 9 Feasible PCP 112.0 0.0 0.0 0.0 100.0 100.0 100.0 10 Feasible OCP 101.0 0.0 0.0 0.0 100.0 100.0 100.0 11 Feasible OCP 101.3 0.0 0.0 0.0 100.0 100.0 100.0 12 Feasible PCP 92.9 0.0 0.0 0.0 100.0 100.0 100.0 13 Feasible PCP 85.0 0.0 0.0 0.0 100.0 100.0 100.0 14 Feasible PCP 82.0 0.0 0.0 0.0 100.0 100.0 100.0 15 Feasible OCP 137.6 0.0 0.0 0.0 100.0 100.0 100.0 16 Feasible OCP 170.4 0.0 0.0 0.0 100.0 100.0 100.0 17 Feasible OCP 291.8 0.0 0.0 0.0 100.0 99.8 99.9 18 Feasible OCP 146.4 0.0 0.0 0.0 100.0 100.0 100.0 19 Feasible OCP 178.3 0.0 0.0 0.0 100.0 100.0 100.0 20 Feasible OCP 207.5 0.0 0.0 0.0 100.0 100.0 1 100.0 166 Table 6.6: Simulation results using the hedging < corridor policy Cases Demand Feasibility Policy Structure Total cost X, x 2 X Pt p 2 p 1 Feasible OCP 104.1 0.0 0.0 0.0 100.0 100.0 100.0 2 Feasible OCP 125.0 0.0 0.0 0.0 100.0 100.0 100.0 3 Feasible OCP 175.4 0.0 0.0 0.0 100.0 100.0 100.0 4 Feasible OCP 106.4 0.0 0.0 0.0 100.0 100.0 100.0 5 Feasible OCP 108.8 0.0 0.0 0.0 100.0 100.0 100.0 6 Feasible OCP 111.1 0.0 0.0 0.0 100.0 100.0 100.0 7 Feasible PCP 110.4 0.0 0.0 0.0 100.0 100.0 100.0 8 Feasible PCP 116.7 0.0 0.0 0.0 100.0 100.0 100.0 9 Feasible PCP 123.0 0.0 0.0 0.0 100.0 100.0 100.0 10 Feasible OCP 104.7 0.0 0.0 0.0 100.0 100.0 100.0 11 Feasible OCP 105.2 0.0 0.0 0.0 100.0 100.0 100.0 12 Feasible PCP 89.0 0.0 0.0 0.0 100.0 100.0 100.0 13 Feasible PCP 75.3 0.0 0.0 0.0 100.0 100.0 100.0 14 Feasible PCP 71.5 0.0 0.0 0.0 100.0 100.0 100.0 15 Feasible OCP 131.2 0.0 0.0 0.0 100.0 100.0 100.0 16 Feasible OCP 184.0 0.0 0.0 0.0 100.0 100.0 100.0 17 Feasible OCP 361.4 0.0 0.0 0.0 100.0 99.8 99.9 18 Feasible OCP 160.0 0.0 0.0 0.0 100.0 100.0 100.0 19 Feasible OCP 202.7 0.0 0.0 0.0 100.0 99.9 100.0 20 Feasible OCP 255.5 0.0 0.0 0.0 100.0 100.0 100.0 167 6.7 Summary In this chapter, we studied an unreliable one-machine two-part-type manufacturing system with random nonresumable setups. We formulated the problem as a stochastic optimal control problem. To obtain the optimal setup and production control policy, we used a simple technique due to Kushner and Dupuis (1992) and converted the problem to a classical infinite horizon discrete time discrete state Markov Decision Process. The latter was solved numerically using the value iteration algorithm of dynamic programming. The optimal policy was shown to have two different structures (Orthogonal Corridor Policy or Parallel Corridor Policy) . Based on results of the corresponding deterministic system, a heuristic, implementable in real-time, was proposed. Simulation results show that the heuristic policy performs very well when applied to systems that are sufficiently reliable producing similar part types or not heavily loaded. CHAPTER 7 SUMMARY AND SUGGESTIONS FOR FUTURE RESEARCH 7.1 Summary In this dissertation, a study of a manufacturing system, comprising a single machine that is able to produce different part types or products, has been presented. The manufacturing system is inflexible in the sense that when production is switched from one part type to the next a setup time and setup cost are incurred. The manufacturing system has a finite production capacity, and has to satisfy a constant demand, for each part type, over a certain planning horizon. Two different scenarios has been considered: A manufacturing system with a perfectly reliable machine, which we referred to as the deterministic system, and a manufacturing system with an unreliable machine subject to random failures, random repairs and random setup times, which we referred to as the stochastic system. The study of the deterministic system has been presented in Chapters 2, 3, 4, and 5. In Chapter 2, the problem has been formulated as an optimal control problem. In Chapter 3, the optimal steady state solution of a two-part- type system has been obtained analytically. The problem was then extended to the multi -part-type case. For the latter, a nonlinear programming algorithm, to obtain the optimal steady state solution, was proposed. In Chapters 4 and 5, a study of the transient behavior of the deterministic two-part-type system was presented. In Chapter 4, a system incurring setup times and no setup costs, when production is switched from one part type to the next, was studied. The 168 169 optimal transient solution was obtained by partitioning the production surplus space (inventory /backlog space) of the system into two mutually exclusive major regions. For initial production surplus levels in one region, the optimal transient solution was obtained by inspection. For initial production surplus levels in the other region, the optimal transient solution was first obtained by means of a novel algorithm. Then, based on insights from the numerical solution, the optimal transient solution was derived analytically. In Chapter 5, The results of Chapter 4 were extended to the more general case where both setup times and costs are incurred, when production is switched from one part type to the next. The complete optimal transient solution was obtained analytically as a feedback law. In Chapter 6, a study of the stochastic system was presented. The problem was formulated as a stochastic optimal control problem, and necessary and sufficient optimal conditions were given as a system of nonlinear partial differential equations, known as the Hamilton-Jacobi-Bellman (HJB) optimality partial differential equations. A simple technique, due to Kushner and Dupuis (1992), was then adopted to convert the HJB optimality partial differential equations to a classical infinite horizon discrete time discrete state Markov Decision Process. The optimal production and setup controls were then obtained numerically by means of the classical value iteration algorithm of dynamic programming. The optimal production policy has been shown to be of the bang-bang type, where either a part type is produced at the maximum machine production rate or not produced at all. The optimal setup policy has been shown to have two different structures, depending on the parameters of the system: A Orthogonal Corridor Structure and a Parallel Corridor Structure. Given the computational difficulties of the exact technique, a simple heuristic policy (Hedging Corridor Policy), based on results from the deterministic 170 system, was then proposed to approximate the optimal policy. Simulation results showed that the heuristic policy gives a very good approximation for systems that are sufficiently reliable with similar characteristics of the produced part types. 7.2 Suggestions for Future Work This research is only the start of work in a potentially fruitful arena of real-time production scheduling problems. There are many possible directions that can be taken to enhance and further develop solution methodologies for the complex dynamic setup problem beyond the work of this dissertation. For the deterministic system, a natural extension is to consider a multi- part-type system, we already presented (Chapter 3) a nonlinear programming algorithm to obtain the optimal steady state solution of such a system. To obtain the optimal transient solution, a technique similar to the one developed in this dissertation can be adopted. First, the production surplus space has to be partitionrd into two mutually exclusive major regions. One above the cyclic schedule and one below it. For initial production surplus levels above the cyclic schedule, optimal control techniques can be used to determine the optimal trajectory leading to the cyclic schedule from above. For initial surplus levels below the cyclic schedule, the trajectory leading to the cyclic schedule from below can be obtained by extending the algorithm, developed in Chapter 4, to handle multiple part types. Another possible direction would be to analyze the setup problem for a deterministic system with multiple machines and multiple part types in which it is not possible to devote each machine to the production of a single part type. 171 The stochastic system can be extended to the three-part-type system. The latter is further complicated by the sequencing policy that has to be determined. The same numerical technique, developed in Chapter 6, can be applied to this case. By observing the numerical solution, simple heuristics, implementable in real-time, can be proposed. The Hedging Corridor Policy can also be tested on the three-part-type system. If it shows a good performance, then it can be extended to the multi-part-type system. The advantage of the Hedging Corridor Policy lies in its simplicity. LIST OF REFERENCES Akella, R., Choong, Y. and Gershwin, S. B. (1984), "Performance of Hierarchical Scheduling Policy," IEEE Transaction on Components. Hybrids, and Manufacturing Technology , CHMT-7(3), 225-240. Akella, R. and Kumar, P. R. (1986), "Optimal Control of Production Rate in a Failure Prone Manufacturing System," IEEE Transaction on Automatic Control , AC-3 1(2), 116-126. Arizono, I., Yokoi, S. and Ohta H. (1989), "The Effects of Varying Production Rates on Inventory Control," Journal of Operational Research Society , 40, 789- 796. Axsater, S. (1985), "Performance Bounds for Lot Sizing Heuristics," Management Science . 31, 634-640. Bai, S. X. and Burhanpurwala, J. H. (1993), "Real-Time Scheduling Algorithm for a Dynamic Setup Problem," SPIE'S International Symposium on Optics for Manufacturing and Advanced Automation . Boston. Bazaraa, M. S. and Shetty, C. M. (1979), NonLinear Programming Theory and Algorithms . John Wiley 8s Sons, New York. Bertsekas, D. P. (1987), Dynamic Programming. Deterministic and Stochastic Models . Prentice Hall, Inc., Englewood Cliffs, NJ. Boctor, F. F. (1982), "The Two-Product, Single-Machine, Static-Demand, Infinite Horizon Lot Scheduling Problem," Management Science , 28, 798-807. Bomberger, E. (1966), "A Dynamic Programming Approach to a Lot Scheduling Problem," Management Science . 12, 778-784. Boukas, E. K. and Haurie, A. (1990), "Manufacturing Flow Control and Preventive Maintenance: A Stochastic Control Approach," IEEE transaction on Automatic Control . AC-35(9), 1024-1031. Buzacott, J. A. and Ozkarahan, I. A. (1983), "One- and Two-Stage Scheduling of Two Products with Distributed Inserted Idle Time: The Benefits of a 172 173 Controllable Production Rate," Naval Research Logistics Quarterly , 30, 675- 696. Carreno, J. J. (1990), "Economic Lot Scheduling for Multiple Products on Parallel Identical Processors," Management Science , 36, 348-358. Connolly, S. (1992), "A Real-Time Policy for Performing Setup Changes in a Manufacturing System," Master's Thesis, MIT Operations Research Center, MIT Laboratory for Manufacturing and Productivity. Report LMP-92-0005, Cambridge. Connolly, S., Dallery, Y. and Gershwin, S. B. (1992), "A Real-Time Policy for Performing Setup Changes in a Manufacturing System," Proceedings of the 31st IEEE Conference on Decision and Control , Tucson, Arizona, December 16-18. Dobson, G. (1987), "The Economic Lot Scheduling Problem: Achieving Feasibility Using Time-Varying Lot Sizes," Operations Research , 35, 764-771. Elmaghraby, S. E. (1978), "The Economic Lot Scheduling Problem (ELSP): Review and Extensions," Management Science , 24, 587-598. Gallego, G. (1989), "Real-Time Scheduling of Several Products on a Single Facility With Setups and Backorders," Ph.D. Thesis, School of OR&IE, Cornell University, Ithaca, NY. Gershwin, S. B. (1993), Manufacturing Systems Engineering . Prentice Hall, Inc., Englewood Cliffs, NJ. Goyal, S. K. (1984), "Determination of Economic Production Quantities for a Two-Product Single Machine System," International Journal of Production Research. 22, 121-126. Haurie, A. and L'Ecuyer P. (1986), "Approximation and Bounds in Discrete Event Dynamic Programming," IEEE Transaction on Automatic Control , AC- 31(3), 227-235. Hsu, W. L. (1983), "On the General Feasibility Test of Scheduling Lot Sizes for Several Products on One Machine," Management Science , 29, 93-105. Hu, J. and Caramanis, M. (1995), "Dynamic Set-up Scheduling of Flexible Manufacturing Systems: Design and Stability of Near Optimal General Round Robin Policies," to appear on Discrete Event Systems, IMA Volumes in Mathematics and its Applications Series , Springer- Verlag, New York. Jones, P. and Inman, R. (1989), "When is the Economic Lot Scheduling Problem easy?," HE Transactions , 21, 11-20. 174 Kimemia, J. and Gershwin, S. B. (1983), " An Algorithm for the Computer Control of a flexible Manufacturing System," HE Transactions . 15(4), 353-362. Kushner, H. J. and Dupuis, P. G. (1993), Numerical Methods for Stochastic Control Problems in Continuous Time . Springer- Verlag, New York. Kumar, P. R. and Seidman, T. I. (1990), "Dynamic Instabilities and Stabilization Methods in Distributed Real-Time Scheduling of Manufacturing Systems," IEEE Transactions on Automatic Control . 35, 289-338. Lee, C. Y. and Surya, D. (1989), "Economic Lot Scheduling for Two-Product Problem," HE Transactions, 21, 162-169. Lou, S., Sethi, S. and Sorger, G. (1991), "Analysis of a Class of Real-Time Multiproduct Lot Scheduling Policies," IEEE Transactions on Automatic Control , AC-36(2), 243-248. Lou, S., Sethi, S. and Sorger, G. (1992), "Stability of Real-Time Lot Scheduling Policies for an Unreliable Machine," IEEE Transactions on Automatic Control . AC-37(12), 1966-1970. Maxwell, W. L. and Singh, H. (1983), "The Effect of Restricting Cycle Times in The Economic Lot Scheduling Problem," HE Transactions . 15, 235-241. Moon, I., Gallego, G. and Simchi-Levi, D. (1991), "Controllable Production Rates in a Family Context," International Journal of Production Research . 29, 2459-2470. Newson, E. F. P. (1975a), "Multi-Item Lot Size Scheduling by Heuristic, Part I: With Fixed Resources," Management Science . 21-10, 1186-1193. Newson, E. F. P. (1975b), "Multi-Item Lot Size Scheduling by Heuristic, Part II: With Variable Resources," Management Science . 21-10, 1193-1203. Olsder, G. J. and Suri, R. (1980), "Time Optimal of Parts-Routing in a Manufacturing System with Failure Prone Machines," Proceedings of the 19th IEEE Conference on Decision and Control . Albuquerque, NM, 722-727. Perkins, J. and Kumar, P. R. (1989), "Stable, Distributed, Real-Time Scheduling of Flexible Manufacturing/Assembly/Disassembly Systems," IEEE Transactions on Automatic Control . 34, 139-148. Rishel, R. (1975), "Dynamic Programming and Minimum Principles for Systems with Jump Markov Disturbances," SIAM Journal On Control . 113(2), 338-371. Roundy, R. (1985), "98%-Effective Integer-Ratio Lot Sizing for One-Warehouse Multi-Retailer Systems," Management Science . 31, 1416-1430. 175 Sharifnia, A., Caramanis, M. and Gershwin S. B. (1991), "Dynamic Setup Scheduling and Flow Control in Manufacturing Systems," Discrete Event Dynamic Systems: Theory and Applications . 1, 149-175. Silver, E. (1990), "Deliberately Slowing Down Output in a Family Production Context," International Journal of Production Research . 28, 17-27. Sivazlian, B. D. and Stanfel, L. E. (1975), Analysis of Systems in Operations Research . Prentice Hall, Inc., Englewood Cliffs, NJ. Srivatsan, N. and Gershwin, S. B. (1990), "Selection of Setup Times in a Hierarchically Controlled Manufacturing System," Proceedings of the 29th IEEE Conference on Decision and Control . Honolulu, HI, December. Tsitsiklis, J. N. (1982), "Characterization of Optimal Policies in a Dynamic Routing Problem," M.I.T. Laboratory for Information and Decision Systems Report No. LIDS-r-1 178, Cambridge. Wagner, H. M. and Whitin, T. M. (1958), "Dynamic Version of the Economic Lot Size Model," Management Science . 5(1), 56-63. Zipkin, P. (1991), "Computing Optimal Lot Sizes in the Economic Lot Scheduling Problem," Operations Research . 39, 53-63. BIOGRAPHICAL SKETCH Mohsen El Hafsi was born January 7, 1963, in Tunis, Tunisia. He received primary and secondary education in Tunis, Tunisia. In June of 1982, he graduated from high school and received his baccalaureate diploma, with high honors, from Lycee Bab El Khadra de Tunis. He attended the Ecole Nationale d'Ingenieurs de Tunis, Tunisia, from 1982 to 1988, receiving the Ingenieur Principal's (Qualified Engineer) degree in industrial engineering, with high honors, in July 1988. From October 1988 to August 1990, he worked as a research engineer in the Institut Regional des Sciences Informatiques et des Telecommunications. In August 1990, he was awarded a scholarship from the Tunisian government to pursue his graduate studies in the United States of America. In January 1991, he joined the Industrial and Systems Engineering Department at the University of Florida to pursue a Master of Science and Doctor of Philosophy degrees. 176 I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Boghos D. Sivaz#an, CI Professor of Industrial/ Systems Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Sherman X. Bai, Cochairman Assistant Professor of Industrial and Systems Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Donald J. arzinga ProfessorKn Industrial and Systems Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Chung-Yee Lee Professor of Industrial and Systems Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. ^ t /fr cJ\ /U*r£ _ Mark C. K. Yang Professor of Statis€icfe This dissertation was submitted to the Graduate Faculty of the College of Engineering and to the Graduate School and was accepted as partial fulfillment of the requirements for the degree of Doctor of Philosophy. August, 1995 Winfred M. Phillips Dean, College of Engineering Karen A. Holbrook Dean, Graduate School